From: Martin Kronbichler <kronbichler@lnm.mw.tum.de> Date: Mon, 28 Apr 2014 19:53:24 +0000 (+0000) Subject: Cleanup file fe_trace.cc X-Git-Tag: v8.2.0-rc1~539 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=39c4e010e3a9e03882a81e5d63937e37ca26b628;p=dealii.git Cleanup file fe_trace.cc git-svn-id: https://svn.dealii.org/trunk@32853 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/fe/fe_trace.h b/deal.II/include/deal.II/fe/fe_trace.h index d098aa16b7..86f2be75be 100644 --- a/deal.II/include/deal.II/fe/fe_trace.h +++ b/deal.II/include/deal.II/fe/fe_trace.h @@ -24,7 +24,7 @@ DEAL_II_NAMESPACE_OPEN /** - * A finite element, which is the trace of Fe_Q elements, that is + * 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 @@ -90,6 +90,13 @@ public: virtual bool has_support_on_face (const unsigned int shape_index, const unsigned int face_index) const; + /** + * Returns a list of constant modes of the element. For this element, it + * simply returns one row with all entries set to true. + */ + virtual std::pair<Table<2,bool>, std::vector<unsigned int> > + get_constant_modes () const; + private: /** * Return vector with dofs per diff --git a/deal.II/source/fe/fe_trace.cc b/deal.II/source/fe/fe_trace.cc index 902ae85f0d..af79e48461 100644 --- a/deal.II/source/fe/fe_trace.cc +++ b/deal.II/source/fe/fe_trace.cc @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2000 - 2013 by the deal.II authors +// Copyright (C) 2000 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -14,98 +14,18 @@ // // --------------------------------------------------------------------- -#include <deal.II/grid/tria.h> -#include <deal.II/dofs/dof_handler.h> -#include <deal.II/grid/grid_generator.h> -#include <deal.II/grid/tria_accessor.h> -#include <deal.II/grid/tria_iterator.h> -#include <deal.II/dofs/dof_accessor.h> - #include <deal.II/base/config.h> #include <deal.II/base/tensor_product_polynomials.h> #include <deal.II/fe/fe_poly_face.h> - -#include <deal.II/fe/fe_q.h> #include <deal.II/fe/fe_tools.h> -#include <deal.II/dofs/dof_tools.h> -#include <deal.II/fe/fe_values.h> -#include <deal.II/base/quadrature_lib.h> -#include <deal.II/base/function.h> -#include <deal.II/numerics/vector_tools.h> -#include <deal.II/numerics/matrix_tools.h> -#include <deal.II/lac/vector.h> -#include <deal.II/lac/full_matrix.h> -#include <deal.II/lac/sparse_matrix.h> -#include <deal.II/lac/compressed_sparsity_pattern.h> -#include <deal.II/lac/solver_cg.h> -#include <deal.II/lac/precondition.h> -#include <deal.II/numerics/data_out.h> -#include <deal.II/numerics/data_out_faces.h> -#include <fstream> -#include <iostream> - +#include <deal.II/fe/fe_trace.h> -#include <deal.II/fe/fe_face.h> -#include <deal.II/fe/fe_poly_face.templates.h> #include <sstream> +#include <fstream> +#include <iostream> DEAL_II_NAMESPACE_OPEN -/** FE-Class for continuous face functions -* adaption of FE_Q and FE_FaceQ -*/ - -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); -}; template <int dim, int spacedim> @@ -246,6 +166,21 @@ FE_TraceQ<dim,spacedim>::has_support_on_face (const unsigned int shape_index, return false; } + + +template <int dim, int spacedim> +std::pair<Table<2,bool>, std::vector<unsigned int> > +FE_TraceQ<dim,spacedim>::get_constant_modes () const +{ + Table<2,bool> constant_modes(1, this->dofs_per_cell); + for (unsigned int i=0; i<this->dofs_per_cell; ++i) + constant_modes(0,i) = true; + return std::pair<Table<2,bool>, std::vector<unsigned int> > + (constant_modes, std::vector<unsigned int>(1, 0)); +} + + + template <int dim, int spacedim> std::vector<unsigned int> FE_TraceQ<dim,spacedim>::get_dpo_vector (const unsigned int deg)