17 #ifndef __deal2__fe_system_h
18 #define __deal2__fe_system_h
24 #include <deal.II/base/config.h>
25 #include <deal.II/base/thread_management.h>
26 #include <deal.II/fe/fe.h>
155 template <
int dim,
int spacedim=dim>
174 const unsigned int n_elements);
221 const std::vector<unsigned int> &multiplicities);
236 virtual std::string
get_name ()
const;
263 const unsigned int component)
const;
292 const unsigned int component)
const;
323 const unsigned int component)
const;
358 const unsigned int face_index)
const;
453 const unsigned int face,
454 const bool face_orientation =
true,
455 const bool face_flip =
false,
456 const bool face_rotation =
false)
const;
517 const unsigned int subface,
538 std::vector<std::pair<unsigned int, unsigned int> >
546 std::vector<std::pair<unsigned int, unsigned int> >
554 std::vector<std::pair<unsigned int, unsigned int> >
626 CellSimilarity::Similarity &cell_similarity)
const;
637 const unsigned int face_no,
652 const unsigned int face_no,
653 const unsigned int sub_no,
673 const unsigned int face_no,
674 const unsigned int sub_no,
676 CellSimilarity::Similarity cell_similarity,
702 std::vector<std::pair<std_cxx1x::shared_ptr<const FiniteElement<dim,spacedim> >,
733 const unsigned int N1,
735 const unsigned int N2=0,
737 const unsigned int N3=0,
739 const unsigned int N4=0,
741 const unsigned int N5=0);
748 const std::vector<unsigned int> &multiplicities);
758 static std::vector<bool>
761 const unsigned int N1,
763 const unsigned int N2=0,
765 const unsigned int N3=0,
767 const unsigned int N4=0,
769 const unsigned int N5=0);
776 static std::vector<bool>
779 const std::vector<unsigned int> &multiplicities);
785 static std::vector<ComponentMask>
787 const unsigned int N1,
789 const unsigned int N2=0,
791 const unsigned int N3=0,
793 const unsigned int N4=0,
795 const unsigned int N5=0);
802 static std::vector<ComponentMask>
804 const std::vector<unsigned int> &multiplicities);
812 const std::vector<unsigned int> &multiplicities);
834 template <
int structdim>
835 std::vector<std::pair<unsigned int, unsigned int> >
925 typename std::vector<typename FiniteElement<dim,spacedim>::InternalDataBase *>
base_fe_datas;
945 DEAL_II_NAMESPACE_CLOSE
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
static const unsigned int invalid_unsigned_int
static FiniteElementData< dim > multiply_dof_numbers(const FiniteElement< dim, spacedim > *fe1, const unsigned int N1, const FiniteElement< dim, spacedim > *fe2=NULL, const unsigned int N2=0, const FiniteElement< dim, spacedim > *fe3=NULL, const unsigned int N3=0, const FiniteElement< dim, spacedim > *fe4=NULL, const unsigned int N4=0, const FiniteElement< dim, spacedim > *fe5=NULL, const unsigned int N5=0)
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual std::string get_name() const
void compute_fill(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim_1 > &quadrature, CellSimilarity::Similarity cell_similarity, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename Mapping< dim, spacedim >::InternalDataBase &fe_data, FEValuesData< dim, spacedim > &data) const
virtual bool hp_constraints_are_implemented() const
std::vector< FEValuesData< dim, spacedim > * > base_fe_values_datas
void initialize_unit_face_support_points()
FESystem(const FiniteElement< dim, spacedim > &fe, const unsigned int n_elements)
virtual void fill_fe_subface_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename Mapping< dim, spacedim >::InternalDataBase &fe_data, FEValuesData< dim, spacedim > &data) const
void delete_fe_values_data(const unsigned int base_no)
void set_fe_data(const unsigned int base_no, typename FiniteElement< dim, spacedim >::InternalDataBase *)
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual std::size_t memory_consumption() const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
void initialize(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities)
virtual UpdateFlags update_each(const UpdateFlags flags) const
virtual Point< dim > unit_support_point(const unsigned int index) const
virtual void fill_fe_face_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename Mapping< dim, spacedim >::InternalDataBase &fe_data, FEValuesData< dim, spacedim > &data) const
virtual FiniteElement< dim, spacedim > * clone() const
virtual Mapping< dim, spacedim >::InternalDataBase * get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
void set_fe_values_data(const unsigned int base_no, FEValuesData< dim, spacedim > *)
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
static std::vector< ComponentMask > compute_nonzero_components(const FiniteElement< dim, spacedim > *fe1, const unsigned int N1, const FiniteElement< dim, spacedim > *fe2=NULL, const unsigned int N2=0, const FiniteElement< dim, spacedim > *fe3=NULL, const unsigned int N3=0, const FiniteElement< dim, spacedim > *fe4=NULL, const unsigned int N4=0, const FiniteElement< dim, spacedim > *fe5=NULL, const unsigned int N5=0)
static std::vector< bool > compute_restriction_is_additive_flags(const FiniteElement< dim, spacedim > *fe1, const unsigned int N1, const FiniteElement< dim, spacedim > *fe2=NULL, const unsigned int N2=0, const FiniteElement< dim, spacedim > *fe3=NULL, const unsigned int N3=0, const FiniteElement< dim, spacedim > *fe4=NULL, const unsigned int N4=0, const FiniteElement< dim, spacedim > *fe5=NULL, const unsigned int N5=0)
static const unsigned int invalid_face_number
virtual void clear_first_cell()
void initialize_unit_support_points()
void initialize_quad_dof_index_permutation()
virtual Point< dim-1 > unit_face_support_point(const unsigned int index) const
std::vector< std::pair< std_cxx1x::shared_ptr< const FiniteElement< dim, spacedim > >, unsigned int > > base_elements
std::vector< typename FiniteElement< dim, spacedim >::InternalDataBase * > base_fe_datas
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false) const
virtual void fill_fe_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename Mapping< dim, spacedim >::InternalDataBase &fe_data, FEValuesData< dim, spacedim > &data, CellSimilarity::Similarity &cell_similarity) const
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
FiniteElement< dim, spacedim >::InternalDataBase & get_fe_data(const unsigned int base_no) const
unsigned int n_base_elements() const
InternalData(const unsigned int n_base_elements, const bool compute_hessians)
virtual UpdateFlags update_once(const UpdateFlags flags) const
void build_interface_constraints()
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const
const bool compute_hessians
FEValuesData< dim, spacedim > & get_fe_values_data(const unsigned int base_no) const
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature) const
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
std::vector< std::pair< unsigned int, unsigned int > > hp_object_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const