Reference documentation for deal.II version 8.1.0
fe_system.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_system.h 30333 2013-08-18 03:23:21Z bangerth @f$
3 //
4 // Copyright (C) 1999 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__fe_system_h
18 #define __deal2__fe_system_h
19 
20 
21 /*---------------------------- fe_system.h ---------------------------*/
22 
23 
24 #include <deal.II/base/config.h>
25 #include <deal.II/base/thread_management.h>
26 #include <deal.II/fe/fe.h>
27 #include <vector>
28 #include <utility>
29 
31 
155 template <int dim, int spacedim=dim>
156 class FESystem : public FiniteElement<dim,spacedim>
157 {
158 public:
159 
174  const unsigned int n_elements);
175 
181  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
182  const FiniteElement<dim,spacedim> &fe2, const unsigned int n2);
183 
189  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
190  const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
191  const FiniteElement<dim,spacedim> &fe3, const unsigned int n3);
192 
198  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
199  const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
200  const FiniteElement<dim,spacedim> &fe3, const unsigned int n3,
201  const FiniteElement<dim,spacedim> &fe4, const unsigned int n4);
202 
208  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
209  const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
210  const FiniteElement<dim,spacedim> &fe3, const unsigned int n3,
211  const FiniteElement<dim,spacedim> &fe4, const unsigned int n4,
212  const FiniteElement<dim,spacedim> &fe5, const unsigned int n5);
213 
220  FESystem (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
221  const std::vector<unsigned int> &multiplicities);
222 
226  virtual ~FESystem ();
227 
236  virtual std::string get_name () const;
237 
250  virtual double shape_value (const unsigned int i,
251  const Point<dim> &p) const;
252 
261  virtual double shape_value_component (const unsigned int i,
262  const Point<dim> &p,
263  const unsigned int component) const;
264 
279  virtual Tensor<1,dim> shape_grad (const unsigned int i,
280  const Point<dim> &p) const;
281 
290  virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
291  const Point<dim> &p,
292  const unsigned int component) const;
293 
308  virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
309  const Point<dim> &p) const;
310 
319  virtual
321  shape_grad_grad_component (const unsigned int i,
322  const Point<dim> &p,
323  const unsigned int component) const;
324 
336  virtual void
338  FullMatrix<double> &matrix) const;
339 
346  virtual const FiniteElement<dim,spacedim> &
347  base_element (const unsigned int index) const;
348 
357  virtual bool has_support_on_face (const unsigned int shape_index,
358  const unsigned int face_index) const;
359 
382  virtual const FullMatrix<double> &
383  get_restriction_matrix (const unsigned int child,
385 
412  virtual const FullMatrix<double> &
413  get_prolongation_matrix (const unsigned int child,
415 
451  virtual
452  unsigned int face_to_cell_index (const unsigned int face_dof_index,
453  const unsigned int face,
454  const bool face_orientation = true,
455  const bool face_flip = false,
456  const bool face_rotation = false) const;
457 
461  virtual
462  Point<dim>
463  unit_support_point (const unsigned int index) const;
464 
468  virtual
469  Point<dim-1>
470  unit_face_support_point (const unsigned int index) const;
471 
484  virtual bool hp_constraints_are_implemented () const;
485 
498  virtual void
500  FullMatrix<double> &matrix) const;
501 
502 
515  virtual void
517  const unsigned int subface,
518  FullMatrix<double> &matrix) const;
519 
537  virtual
538  std::vector<std::pair<unsigned int, unsigned int> >
540 
545  virtual
546  std::vector<std::pair<unsigned int, unsigned int> >
547  hp_line_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
548 
553  virtual
554  std::vector<std::pair<unsigned int, unsigned int> >
555  hp_quad_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
556 
565  virtual
569 
578  virtual std::size_t memory_consumption () const;
579 
580 protected:
584  virtual UpdateFlags update_once (const UpdateFlags flags) const;
585 
589  virtual UpdateFlags update_each (const UpdateFlags flags) const;
590 
596  virtual FiniteElement<dim,spacedim> *clone() const;
597 
598  virtual typename Mapping<dim,spacedim>::InternalDataBase *
599  get_data (const UpdateFlags update_flags,
600  const Mapping<dim,spacedim> &mapping,
601  const Quadrature<dim> &quadrature) const ;
602 
603  virtual typename Mapping<dim,spacedim>::InternalDataBase *
604  get_face_data (const UpdateFlags update_flags,
605  const Mapping<dim,spacedim> &mapping,
606  const Quadrature<dim-1> &quadrature) const ;
607 
608  virtual typename Mapping<dim,spacedim>::InternalDataBase *
609  get_subface_data (const UpdateFlags update_flags,
610  const Mapping<dim,spacedim> &mapping,
611  const Quadrature<dim-1> &quadrature) const ;
612 
619  virtual void
620  fill_fe_values (const Mapping<dim,spacedim> &mapping,
621  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
622  const Quadrature<dim> &quadrature,
623  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
624  typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
626  CellSimilarity::Similarity &cell_similarity) const;
627 
634  virtual void
636  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
637  const unsigned int face_no,
638  const Quadrature<dim-1> &quadrature,
639  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
640  typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
641  FEValuesData<dim,spacedim> &data) const ;
642 
649  virtual void
651  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
652  const unsigned int face_no,
653  const unsigned int sub_no,
654  const Quadrature<dim-1> &quadrature,
655  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
656  typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
657  FEValuesData<dim,spacedim> &data) const ;
658 
659 
670  template <int dim_1>
671  void compute_fill (const Mapping<dim,spacedim> &mapping,
672  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
673  const unsigned int face_no,
674  const unsigned int sub_no,
675  const Quadrature<dim_1> &quadrature,
676  CellSimilarity::Similarity cell_similarity,
677  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
678  typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
679  FEValuesData<dim,spacedim> &data) const ;
680 
681 private:
682 
687 
702  std::vector<std::pair<std_cxx1x::shared_ptr<const FiniteElement<dim,spacedim> >,
703  unsigned int> >
705 
706 
712 
718 
724 
733  const unsigned int N1,
734  const FiniteElement<dim,spacedim> *fe2=NULL,
735  const unsigned int N2=0,
736  const FiniteElement<dim,spacedim> *fe3=NULL,
737  const unsigned int N3=0,
738  const FiniteElement<dim,spacedim> *fe4=NULL,
739  const unsigned int N4=0,
740  const FiniteElement<dim,spacedim> *fe5=NULL,
741  const unsigned int N5=0);
742 
747  multiply_dof_numbers (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
748  const std::vector<unsigned int> &multiplicities);
749 
750 
751 
758  static std::vector<bool>
760  const FiniteElement<dim,spacedim> *fe1,
761  const unsigned int N1,
762  const FiniteElement<dim,spacedim> *fe2=NULL,
763  const unsigned int N2=0,
764  const FiniteElement<dim,spacedim> *fe3=NULL,
765  const unsigned int N3=0,
766  const FiniteElement<dim,spacedim> *fe4=NULL,
767  const unsigned int N4=0,
768  const FiniteElement<dim,spacedim> *fe5=NULL,
769  const unsigned int N5=0);
770 
776  static std::vector<bool>
778  const std::vector<const FiniteElement<dim,spacedim>*> &fes,
779  const std::vector<unsigned int> &multiplicities);
780 
781 
785  static std::vector<ComponentMask>
787  const unsigned int N1,
788  const FiniteElement<dim,spacedim> *fe2=NULL,
789  const unsigned int N2=0,
790  const FiniteElement<dim,spacedim> *fe3=NULL,
791  const unsigned int N3=0,
792  const FiniteElement<dim,spacedim> *fe4=NULL,
793  const unsigned int N4=0,
794  const FiniteElement<dim,spacedim> *fe5=NULL,
795  const unsigned int N5=0);
796 
802  static std::vector<ComponentMask>
803  compute_nonzero_components (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
804  const std::vector<unsigned int> &multiplicities);
805 
811  void initialize (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
812  const std::vector<unsigned int> &multiplicities);
813 
817  void build_cell_tables();
818 
822  void build_face_tables();
823 
828 
834  template <int structdim>
835  std::vector<std::pair<unsigned int, unsigned int> >
837 
844  class InternalData : public FiniteElement<dim,spacedim>::InternalDataBase
845  {
846  public:
852  InternalData (const unsigned int n_base_elements,
853  const bool compute_hessians);
854 
859  ~InternalData();
860 
864  const bool compute_hessians;
865 
870  void set_fe_data(const unsigned int base_no,
872 
878  get_fe_data (const unsigned int base_no) const;
879 
880 
885  void set_fe_values_data (const unsigned int base_no,
887 
892  FEValuesData<dim,spacedim> &get_fe_values_data (const unsigned int base_no) const;
893 
902  void delete_fe_values_data (const unsigned int base_no);
903 
911  virtual void clear_first_cell ();
912 
913  private:
914 
925  typename std::vector<typename FiniteElement<dim,spacedim>::InternalDataBase *> base_fe_datas;
926 
935  std::vector<FEValuesData<dim,spacedim> *> base_fe_values_datas;
936  };
937 
938  /*
939  * Mutex for protecting initialization of restriction and embedding matrix.
940  */
941  mutable Threads::Mutex mutex;
942 };
943 
944 
945 DEAL_II_NAMESPACE_CLOSE
946 
947 /*---------------------------- fe_system.h ---------------------------*/
948 #endif
949 /*---------------------------- fe_system.h ---------------------------*/
void initialize_unit_support_points()
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)
static const unsigned int invalid_unsigned_int
Definition: types.h:191
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) 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 Tensor< 1, dim > shape_grad(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)
virtual std::string get_name() const
virtual Tensor< 1, dim > shape_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
std::vector< FEValuesData< dim, spacedim > * > base_fe_values_datas
Definition: fe_system.h:935
std::vector< std::pair< unsigned int, unsigned int > > hp_object_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
void delete_fe_values_data(const unsigned int base_no)
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
void initialize_quad_dof_index_permutation()
void set_fe_data(const unsigned int base_no, typename FiniteElement< dim, spacedim >::InternalDataBase *)
virtual Point< dim > unit_support_point(const unsigned int index) const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) 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 double shape_value(const unsigned int i, const Point< dim > &p) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual bool hp_constraints_are_implemented() const
void build_face_tables()
virtual Point< dim-1 > unit_face_support_point(const unsigned int index) const
virtual std::size_t memory_consumption() const
UpdateFlags
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
static const unsigned int invalid_face_number
Definition: fe_system.h:686
virtual FiniteElement< dim, spacedim > * clone() const
void initialize_unit_face_support_points()
virtual ~FESystem()
void set_fe_values_data(const unsigned int base_no, FEValuesData< dim, spacedim > *)
unsigned int n_base_elements() const
Definition: fe.h:2062
virtual UpdateFlags update_once(const UpdateFlags flags) 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 void clear_first_cell()
void initialize(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities)
std::vector< std::pair< std_cxx1x::shared_ptr< const FiniteElement< dim, spacedim > >, unsigned int > > base_elements
Definition: fe_system.h:704
virtual Mapping< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const
std::vector< typename FiniteElement< dim, spacedim >::InternalDataBase * > base_fe_datas
Definition: fe_system.h:925
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) 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 bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
FiniteElement< dim, spacedim >::InternalDataBase & get_fe_data(const unsigned int base_no) const
InternalData(const unsigned int n_base_elements, const bool compute_hessians)
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
FESystem(const FiniteElement< dim, spacedim > &fe, const unsigned int n_elements)
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) 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_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
void build_cell_tables()
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)
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
const bool compute_hessians
Definition: fe_system.h:864
FEValuesData< dim, spacedim > & get_fe_values_data(const unsigned int base_no) const
virtual UpdateFlags update_each(const UpdateFlags flags) const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const