Reference documentation for deal.II version 8.4.1
fe.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2016 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii__fe_h
17 #define dealii__fe_h
18 
19 #include <deal.II/base/config.h>
20 #include <deal.II/base/geometry_info.h>
21 #include <deal.II/fe/fe_base.h>
22 #include <deal.II/fe/fe_values_extractors.h>
23 #include <deal.II/fe/fe_update_flags.h>
24 #include <deal.II/fe/component_mask.h>
25 #include <deal.II/fe/block_mask.h>
26 #include <deal.II/fe/mapping.h>
27 
28 DEAL_II_NAMESPACE_OPEN
29 
30 template <int dim, int spacedim> class FEValuesBase;
31 template <int dim, int spacedim> class FEValues;
32 template <int dim, int spacedim> class FEFaceValues;
33 template <int dim, int spacedim> class FESubfaceValues;
34 template <int dim, int spacedim> class FESystem;
35 namespace hp
36 {
37  template <int dim, int spacedim> class FECollection;
38 }
39 
40 
574 template <int dim, int spacedim=dim>
575 class FiniteElement : public Subscriptor,
576  public FiniteElementData<dim>
577 {
578 public:
582  static const unsigned int space_dimension = spacedim;
583 
611  {
612  private:
617 
618  public:
623  InternalDataBase ();
624 
628  virtual ~InternalDataBase ();
629 
645 
649  virtual std::size_t memory_consumption () const;
650  };
651 
652 public:
695  FiniteElement (const FiniteElementData<dim> &fe_data,
696  const std::vector<bool> &restriction_is_additive_flags,
697  const std::vector<ComponentMask> &nonzero_components);
698 
703  virtual ~FiniteElement ();
704 
711  virtual FiniteElement<dim,spacedim> *clone() const = 0;
712 
723  virtual std::string get_name () const = 0;
724 
746  const FiniteElement<dim,spacedim> &operator[] (const unsigned int fe_index) const;
747 
773  virtual double shape_value (const unsigned int i,
774  const Point<dim> &p) const;
775 
782  virtual double shape_value_component (const unsigned int i,
783  const Point<dim> &p,
784  const unsigned int component) const;
785 
807  virtual Tensor<1,dim> shape_grad (const unsigned int i,
808  const Point<dim> &p) const;
809 
816  virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
817  const Point<dim> &p,
818  const unsigned int component) const;
819 
841  virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
842  const Point<dim> &p) const;
843 
850  virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i,
851  const Point<dim> &p,
852  const unsigned int component) const;
853 
875  virtual Tensor<3,dim> shape_3rd_derivative (const unsigned int i,
876  const Point<dim> &p) const;
877 
884  virtual Tensor<3,dim> shape_3rd_derivative_component (const unsigned int i,
885  const Point<dim> &p,
886  const unsigned int component) const;
887 
909  virtual Tensor<4,dim> shape_4th_derivative (const unsigned int i,
910  const Point<dim> &p) const;
911 
918  virtual Tensor<4,dim> shape_4th_derivative_component (const unsigned int i,
919  const Point<dim> &p,
920  const unsigned int component) const;
931  virtual bool has_support_on_face (const unsigned int shape_index,
932  const unsigned int face_index) const;
933 
935 
955  virtual const FullMatrix<double> &
956  get_restriction_matrix (const unsigned int child,
958 
988  virtual const FullMatrix<double> &
989  get_prolongation_matrix (const unsigned int child,
991 
1013  bool prolongation_is_implemented () const;
1014 
1031 
1053  bool restriction_is_implemented () const;
1054 
1071 
1072 
1081  bool restriction_is_additive (const unsigned int index) const;
1082 
1094  const FullMatrix<double> &constraints (const ::internal::SubfaceCase<dim> &subface_case=::internal::SubfaceCase<dim>::case_isotropic) const;
1095 
1111  bool constraints_are_implemented (const ::internal::SubfaceCase<dim> &subface_case=::internal::SubfaceCase<dim>::case_isotropic) const;
1112 
1113 
1135  virtual bool hp_constraints_are_implemented () const;
1136 
1137 
1149  virtual void
1151  FullMatrix<double> &matrix) const;
1153 
1171  virtual void
1173  FullMatrix<double> &matrix) const;
1174 
1175 
1187  virtual void
1189  const unsigned int subface,
1190  FullMatrix<double> &matrix) const;
1192 
1193 
1214  virtual
1215  std::vector<std::pair<unsigned int, unsigned int> >
1216  hp_vertex_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
1217 
1222  virtual
1223  std::vector<std::pair<unsigned int, unsigned int> >
1224  hp_line_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
1225 
1230  virtual
1231  std::vector<std::pair<unsigned int, unsigned int> >
1232  hp_quad_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
1233 
1243  virtual
1246 
1248 
1258  bool operator == (const FiniteElement<dim,spacedim> &) const;
1259 
1292  std::pair<unsigned int, unsigned int>
1293  system_to_component_index (const unsigned int index) const;
1294 
1304  unsigned int component_to_system_index(const unsigned int component,
1305  const unsigned int index) const;
1306 
1316  std::pair<unsigned int, unsigned int>
1317  face_system_to_component_index (const unsigned int index) const;
1318 
1327  unsigned int adjust_quad_dof_index_for_face_orientation (const unsigned int index,
1328  const bool face_orientation,
1329  const bool face_flip,
1330  const bool face_rotation) const;
1331 
1386  virtual
1387  unsigned int face_to_cell_index (const unsigned int face_dof_index,
1388  const unsigned int face,
1389  const bool face_orientation = true,
1390  const bool face_flip = false,
1391  const bool face_rotation = false) const;
1392 
1400  unsigned int adjust_line_dof_index_for_line_orientation (const unsigned int index,
1401  const bool line_orientation) const;
1402 
1419  const ComponentMask &
1420  get_nonzero_components (const unsigned int i) const;
1421 
1432  unsigned int
1433  n_nonzero_components (const unsigned int i) const;
1434 
1445  bool
1446  is_primitive (const unsigned int i) const;
1447 
1453 
1465  unsigned int n_base_elements () const;
1466 
1471  virtual
1473  base_element (const unsigned int index) const;
1474 
1481  unsigned int
1482  element_multiplicity (const unsigned int index) const;
1483 
1503  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1504  system_to_base_index (const unsigned int index) const;
1505 
1514  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1515  face_system_to_base_index (const unsigned int index) const;
1516 
1521  types::global_dof_index first_block_of_base (const unsigned int b) const;
1522 
1535  std::pair<unsigned int, unsigned int>
1536  component_to_base_index (const unsigned int component) const;
1537 
1538 
1543  std::pair<unsigned int,unsigned int>
1544  block_to_base_index (const unsigned int block) const;
1545 
1549  std::pair<unsigned int,types::global_dof_index>
1550  system_to_block_index (const unsigned int component) const;
1551 
1555  unsigned int
1556  component_to_block_index (const unsigned int component) const;
1557 
1559 
1578  component_mask (const FEValuesExtractors::Scalar &scalar) const;
1579 
1593  component_mask (const FEValuesExtractors::Vector &vector) const;
1594 
1609  component_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1610 
1626  component_mask (const BlockMask &block_mask) const;
1627 
1648  BlockMask
1649  block_mask (const FEValuesExtractors::Scalar &scalar) const;
1650 
1667  BlockMask
1668  block_mask (const FEValuesExtractors::Vector &vector) const;
1669 
1687  BlockMask
1688  block_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1689 
1712  BlockMask
1713  block_mask (const ComponentMask &component_mask) const;
1714 
1730  virtual std::pair<Table<2,bool>,std::vector<unsigned int> >
1731  get_constant_modes () const;
1732 
1734 
1770  const std::vector<Point<dim> > &
1771  get_unit_support_points () const;
1772 
1790  bool has_support_points () const;
1791 
1804  virtual
1805  Point<dim>
1806  unit_support_point (const unsigned int index) const;
1807 
1834  const std::vector<Point<dim-1> > &
1836 
1845  bool has_face_support_points () const;
1846 
1851  virtual
1852  Point<dim-1>
1853  unit_face_support_point (const unsigned int index) const;
1854 
1862  const std::vector<Point<dim> > &
1864 
1873  bool has_generalized_support_points () const;
1874 
1878  const std::vector<Point<dim-1> > &
1879  get_generalized_face_support_points () const;
1880 
1889  bool
1891 
1936  get_associated_geometry_primitive (const unsigned int cell_dof_index) const;
1937 
1947  virtual
1948  void
1949  interpolate(std::vector<double> &local_dofs,
1950  const std::vector<double> &values) const;
1951 
1961  virtual
1962  void
1963  interpolate(std::vector<double> &local_dofs,
1964  const std::vector<Vector<double> > &values,
1965  unsigned int offset = 0) const;
1966 
1971  virtual
1972  void
1973  interpolate(std::vector<double> &local_dofs,
1974  const VectorSlice<const std::vector<std::vector<double> > > &values) const;
1975 
1977 
1986  virtual std::size_t memory_consumption () const;
1987 
1993  DeclException1 (ExcShapeFunctionNotPrimitive,
1994  int,
1995  << "The shape function with index " << arg1
1996  << " is not primitive, i.e. it is vector-valued and "
1997  << "has more than one non-zero vector component. This "
1998  << "function cannot be called for these shape functions. "
1999  << "Maybe you want to use the same function with the "
2000  << "_component suffix?");
2006  DeclException0 (ExcFENotPrimitive);
2012  DeclExceptionMsg (ExcUnitShapeValuesDoNotExist,
2013  "You are trying to access the values or derivatives of shape functions "
2014  "on the reference cell of an element that does not define its shape "
2015  "functions through mapping from the reference cell. Consequently, "
2016  "you cannot ask for shape function values or derivatives on the "
2017  "reference cell.");
2018 
2025  DeclExceptionMsg (ExcFEHasNoSupportPoints,
2026  "You are trying to access the support points of a finite "
2027  "element that either has no support points at all, or for "
2028  "which the corresponding tables have not been implemented.");
2029 
2036  DeclExceptionMsg (ExcEmbeddingVoid,
2037  "You are trying to access the matrices that describe how "
2038  "to embed a finite element function on one cell into the "
2039  "finite element space on one of its children (i.e., the "
2040  "'embedding' or 'prolongation' matrices). However, the "
2041  "current finite element can either not define this sort of "
2042  "operation, or it has not yet been implemented.");
2043 
2051  DeclExceptionMsg (ExcProjectionVoid,
2052  "You are trying to access the matrices that describe how "
2053  "to restrict a finite element function from the children "
2054  "of one cell to the finite element space defined on their "
2055  "parent (i.e., the 'restriction' or 'projection' matrices). "
2056  "However, the current finite element can either not define "
2057  "this sort of operation, or it has not yet been "
2058  "implemented.");
2059 
2064  DeclException2 (ExcWrongInterfaceMatrixSize,
2065  int, int,
2066  << "The interface matrix has a size of " << arg1
2067  << "x" << arg2
2068  << ", which is not reasonable for the current element "
2069  "in the present dimension.");
2074  DeclException0 (ExcInterpolationNotImplemented);
2075 
2076 protected:
2077 
2091  void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false,
2092  const bool isotropic_prolongation_only=false);
2093 
2106  std::vector<std::vector<FullMatrix<double> > > restriction;
2107 
2120  std::vector<std::vector<FullMatrix<double> > > prolongation;
2121 
2133 
2144  std::vector<Point<dim> > unit_support_points;
2145 
2151  std::vector<Point<dim-1> > unit_face_support_points;
2152 
2157  std::vector<Point<dim> > generalized_support_points;
2158 
2164 
2181 
2196 
2200  std::vector<std::pair<unsigned int, unsigned int> > system_to_component_table;
2201 
2211  std::vector<std::pair<unsigned int, unsigned int> > face_system_to_component_table;
2212 
2229  std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
2231 
2235  std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
2237 
2243 
2259  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int> >
2261 
2267  const std::vector<bool> restriction_is_additive_flags;
2268 
2276  const std::vector<ComponentMask> nonzero_components;
2277 
2285  const std::vector<unsigned int> n_nonzero_components_table;
2286 
2299  interface_constraints_size () const;
2300 
2306  static
2307  std::vector<unsigned int>
2308  compute_n_nonzero_components (const std::vector<ComponentMask> &nonzero_components);
2309 
2330  virtual
2331  UpdateFlags
2332  requires_update_flags (const UpdateFlags update_flags) const = 0;
2333 
2410  virtual
2411  InternalDataBase *
2412  get_data (const UpdateFlags update_flags,
2413  const Mapping<dim,spacedim> &mapping,
2414  const Quadrature<dim> &quadrature,
2416 
2458  virtual
2459  InternalDataBase *
2460  get_face_data (const UpdateFlags update_flags,
2461  const Mapping<dim,spacedim> &mapping,
2462  const Quadrature<dim-1> &quadrature,
2464 
2506  virtual
2507  InternalDataBase *
2508  get_subface_data (const UpdateFlags update_flags,
2509  const Mapping<dim,spacedim> &mapping,
2510  const Quadrature<dim-1> &quadrature,
2512 
2593  virtual
2594  void
2596  const CellSimilarity::Similarity cell_similarity,
2597  const Quadrature<dim> &quadrature,
2598  const Mapping<dim,spacedim> &mapping,
2599  const typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
2600  const ::internal::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
2601  const InternalDataBase &fe_internal,
2603 
2646  virtual
2647  void
2649  const unsigned int face_no,
2650  const Quadrature<dim-1> &quadrature,
2651  const Mapping<dim,spacedim> &mapping,
2652  const typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
2653  const ::internal::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
2654  const InternalDataBase &fe_internal,
2656 
2702  virtual
2703  void
2705  const unsigned int face_no,
2706  const unsigned int sub_no,
2707  const Quadrature<dim-1> &quadrature,
2708  const Mapping<dim,spacedim> &mapping,
2709  const typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
2710  const ::internal::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
2711  const InternalDataBase &fe_internal,
2713 
2714  friend class InternalDataBase;
2715  friend class FEValuesBase<dim,spacedim>;
2716  friend class FEValues<dim,spacedim>;
2717  friend class FEFaceValues<dim,spacedim>;
2718  friend class FESubfaceValues<dim,spacedim>;
2719  friend class FESystem<dim,spacedim>;
2720 };
2721 
2722 
2723 //----------------------------------------------------------------------//
2724 
2725 
2726 template <int dim, int spacedim>
2727 inline
2729 FiniteElement<dim,spacedim>::operator[] (const unsigned int fe_index) const
2730 {
2731  (void)fe_index;
2732  Assert (fe_index == 0,
2733  ExcMessage ("A fe_index of zero is the only index allowed here"));
2734  return *this;
2735 }
2736 
2737 
2738 
2739 template <int dim, int spacedim>
2740 inline
2741 std::pair<unsigned int,unsigned int>
2743 {
2744  Assert (index < system_to_component_table.size(),
2745  ExcIndexRange(index, 0, system_to_component_table.size()));
2746  Assert (is_primitive (index),
2748  return system_to_component_table[index];
2749 }
2750 
2751 
2752 
2753 template <int dim, int spacedim>
2754 inline
2755 unsigned int
2757 {
2758  return base_to_block_indices.size();
2759 }
2760 
2761 
2762 
2763 template <int dim, int spacedim>
2764 inline
2765 unsigned int
2767 {
2768  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
2769 }
2770 
2771 
2772 
2773 template <int dim, int spacedim>
2774 inline
2775 unsigned int
2777  const unsigned int index) const
2778 {
2779  AssertIndexRange(component, this->n_components());
2780  const std::vector<std::pair<unsigned int, unsigned int> >::const_iterator
2781  it = std::find(system_to_component_table.begin(), system_to_component_table.end(),
2782  std::pair<unsigned int, unsigned int>(component, index));
2783 
2784  Assert(it != system_to_component_table.end(),
2785  ExcMessage ("You are asking for the number of the shape function "
2786  "within a system element that corresponds to vector "
2787  "component " + Utilities::int_to_string(component) + " and within this to "
2788  "index " + Utilities::int_to_string(index) + ". But no such "
2789  "shape function exists."));
2790  return std::distance(system_to_component_table.begin(), it);
2791 }
2792 
2793 
2794 
2795 template <int dim, int spacedim>
2796 inline
2797 std::pair<unsigned int,unsigned int>
2799 {
2800  Assert(index < face_system_to_component_table.size(),
2801  ExcIndexRange(index, 0, face_system_to_component_table.size()));
2802 
2803  // in debug mode, check whether the
2804  // function is primitive, since
2805  // otherwise the result may have no
2806  // meaning
2807  //
2808  // since the primitivity tables are
2809  // all geared towards cell dof
2810  // indices, rather than face dof
2811  // indices, we have to work a
2812  // little bit...
2813  //
2814  // in 1d, the face index is equal
2815  // to the cell index
2816  Assert (is_primitive(this->face_to_cell_index(index, 0)),
2818 
2819  return face_system_to_component_table[index];
2820 }
2821 
2822 
2823 
2824 template <int dim, int spacedim>
2825 inline
2826 std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
2828 {
2829  Assert (index < system_to_base_table.size(),
2830  ExcIndexRange(index, 0, system_to_base_table.size()));
2831  return system_to_base_table[index];
2832 }
2833 
2834 
2835 
2836 
2837 template <int dim, int spacedim>
2838 inline
2839 std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
2841 {
2842  Assert(index < face_system_to_base_table.size(),
2843  ExcIndexRange(index, 0, face_system_to_base_table.size()));
2844  return face_system_to_base_table[index];
2845 }
2846 
2847 
2848 
2849 template <int dim, int spacedim>
2850 inline
2853 {
2854  return base_to_block_indices.block_start(index);
2855 }
2856 
2857 
2858 
2859 template <int dim, int spacedim>
2860 inline
2861 std::pair<unsigned int,unsigned int>
2863 {
2864  Assert(index < component_to_base_table.size(),
2865  ExcIndexRange(index, 0, component_to_base_table.size()));
2866 
2867  return component_to_base_table[index].first;
2868 }
2869 
2870 
2871 
2872 template <int dim, int spacedim>
2873 inline
2874 std::pair<unsigned int,unsigned int>
2876 {
2877  return base_to_block_indices.global_to_local(index);
2878 }
2879 
2880 
2881 
2882 template <int dim, int spacedim>
2883 inline
2884 std::pair<unsigned int,types::global_dof_index>
2886 {
2887  Assert (index < this->dofs_per_cell,
2888  ExcIndexRange(index, 0, this->dofs_per_cell));
2889  // The block is computed simply as
2890  // first block of this base plus
2891  // the index within the base blocks
2892  return std::pair<unsigned int, types::global_dof_index>(
2893  first_block_of_base(system_to_base_table[index].first.first)
2894  + system_to_base_table[index].first.second,
2895  system_to_base_table[index].second);
2896 }
2897 
2898 
2899 
2900 template <int dim, int spacedim>
2901 inline
2902 bool
2904 {
2905  Assert(index < this->dofs_per_cell,
2906  ExcIndexRange(index, 0, this->dofs_per_cell));
2907  return restriction_is_additive_flags[index];
2908 }
2909 
2910 
2911 
2912 template <int dim, int spacedim>
2913 inline
2914 const ComponentMask &
2916 {
2917  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
2918  return nonzero_components[i];
2919 }
2920 
2921 
2922 
2923 template <int dim, int spacedim>
2924 inline
2925 unsigned int
2927 {
2928  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
2929  return n_nonzero_components_table[i];
2930 }
2931 
2932 
2933 
2934 template <int dim, int spacedim>
2935 inline
2936 bool
2937 FiniteElement<dim,spacedim>::is_primitive (const unsigned int i) const
2938 {
2939  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
2940 
2941  // return primitivity of a shape
2942  // function by checking whether it
2943  // has more than one non-zero
2944  // component or not. we could cache
2945  // this value in an array of bools,
2946  // but accessing a bit-vector (as
2947  // std::vector<bool> is) is
2948  // probably more expensive than
2949  // just comparing against 1
2950  //
2951  // for good measure, short circuit the test
2952  // if the entire FE is primitive
2953  return (is_primitive() ||
2954  (n_nonzero_components_table[i] == 1));
2955 }
2956 
2957 
2958 
2959 template <int dim, int spacedim>
2960 inline
2963 {
2964  Assert (cell_dof_index < this->dofs_per_cell,
2965  ExcIndexRange (cell_dof_index, 0, this->dofs_per_cell));
2966 
2967  // just go through the usual cases, taking into account how DoFs
2968  // are enumerated on the reference cell
2969  if (cell_dof_index < this->first_line_index)
2970  return GeometryPrimitive::vertex;
2971  else if (cell_dof_index < this->first_quad_index)
2972  return GeometryPrimitive::line;
2973  else if (cell_dof_index < this->first_hex_index)
2974  return GeometryPrimitive::quad;
2975  else
2976  return GeometryPrimitive::hex;
2977 }
2978 
2979 
2980 
2981 DEAL_II_NAMESPACE_CLOSE
2982 
2983 #endif
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:164
const unsigned int first_hex_index
Definition: fe_base.h:259
bool prolongation_is_implemented() const
Definition: fe.cc:673
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual void interpolate(std::vector< double > &local_dofs, const std::vector< double > &values) const
Definition: fe.cc:1103
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const =0
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
Definition: fe.cc:1091
const std::vector< Point< dim > > & get_unit_support_points() const
Definition: fe.cc:956
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2106
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
Definition: fe.h:2885
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:913
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2157
FullMatrix< double > interface_constraints
Definition: fe.h:2132
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:61
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
Definition: fe.h:2776
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:233
bool isotropic_restriction_is_implemented() const
Definition: fe.cc:751
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:851
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:777
::ExceptionBase & ExcMessage(std::string arg1)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1081
Auxiliary class aiding in the handling of block structures like in BlockVector or FESystem...
Definition: block_indices.h:54
Table< 2, int > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2180
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:924
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2260
bool is_primitive() const
virtual bool hp_constraints_are_implemented() const
Definition: fe.cc:789
bool isotropic_prolongation_is_implemented() const
Definition: fe.cc:725
TableIndices< 2 > interface_constraints_size() const
Definition: fe.cc:824
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:210
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:1193
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > face_system_to_base_table
Definition: fe.h:2236
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
Definition: fe_system.cc:1298
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
Definition: fe.cc:1080
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2285
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:443
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:902
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2144
bool restriction_is_implemented() const
Definition: fe.cc:699
unsigned int element_multiplicity(const unsigned int index) const
Definition: fe.h:2766
size_type block_start(const unsigned int i) const
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
Definition: fe.cc:326
const unsigned int first_quad_index
Definition: fe_base.h:254
virtual ~FiniteElement()
Definition: fe.cc:156
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2120
unsigned int global_dof_index
Definition: types.h:88
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
Definition: fe.h:2875
DeclException1(ExcShapeFunctionNotPrimitive, int,<< "The shape function with index "<< arg1<< " is not primitive, i.e. it is vector-valued and "<< "has more than one non-zero vector component. This "<< "function cannot be called for these shape functions. "<< "Maybe you want to use the same function with the "<< "_component suffix?")
#define Assert(cond, exc)
Definition: exceptions.h:294
unsigned int component_to_block_index(const unsigned int component) const
Definition: fe.cc:348
DeclException0(ExcFENotPrimitive)
UpdateFlags
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:198
unsigned int n_nonzero_components(const unsigned int i) const
Definition: fe.h:2926
bool operator==(const FiniteElement< dim, spacedim > &) const
Definition: fe.cc:945
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:175
unsigned int adjust_quad_dof_index_for_face_orientation(const unsigned int index, const bool face_orientation, const bool face_flip, const bool face_rotation) const
Definition: fe.cc:618
Abstract base class for mapping classes.
Definition: dof_tools.h:52
std::vector< Point< dim-1 > > unit_face_support_points
Definition: fe.h:2151
Definition: fe.h:34
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
Definition: fe.h:2862
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
DeclExceptionMsg(ExcUnitShapeValuesDoNotExist,"You are trying to access the values or derivatives of shape functions ""on the reference cell of an element that does not define its shape ""functions through mapping from the reference cell. Consequently, ""you cannot ask for shape function values or derivatives on the ""reference cell.")
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
Definition: fe.cc:650
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
Definition: fe.h:2962
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:2915
virtual Point< dim > unit_support_point(const unsigned int index) const
Definition: fe.cc:1005
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:798
virtual std::string get_name() const =0
unsigned int n_base_elements() const
Definition: fe.h:2756
bool has_generalized_support_points() const
Definition: fe.cc:996
virtual Point< dim-1 > unit_face_support_point(const unsigned int index) const
Definition: fe.cc:1067
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index) const
Definition: fe.h:2798
virtual std::size_t memory_consumption() const
Definition: fe.cc:1173
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
Definition: fe.h:2827
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:187
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:82
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:868
virtual ~InternalDataBase()
Definition: fe.cc:45
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
Definition: fe.cc:885
Definition: hp.h:102
size_type block_size(const unsigned int i) const
const FiniteElement< dim, spacedim > & operator[](const unsigned int fe_index) const
Definition: fe.h:2729
const unsigned int dofs_per_cell
Definition: fe_base.h:283
std::vector< Point< dim-1 > > generalized_face_support_points
Definition: fe.h:2163
DeclException2(ExcWrongInterfaceMatrixSize, int, int,<< "The interface matrix has a size of "<< arg1<< "x"<< arg2<< ", which is not reasonable for the current element ""in the present dimension.")
const std::vector< Point< dim > > & get_generalized_support_points() const
Definition: fe.cc:981
bool restriction_is_additive(const unsigned int index) const
Definition: fe.h:2903
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:244
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
Definition: fe.cc:529
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
Definition: fe.cc:306
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > face_system_to_base_index(const unsigned int index) const
Definition: fe.h:2840
virtual std::size_t memory_consumption() const
Definition: fe.cc:52
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:2742
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValues::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:935
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2200
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:256
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:361
Definition: fe.h:31
const unsigned int first_line_index
Definition: fe_base.h:249
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2276
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Definition: fe.cc:278
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:267
bool has_generalized_face_support_points() const
Definition: fe.cc:1058
unsigned int size() const
virtual InternalDataBase * get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
Definition: fe.cc:1208
std::vector< std::pair< unsigned int, unsigned int > > face_system_to_component_table
Definition: fe.h:2211
const std::vector< Point< dim-1 > > & get_unit_face_support_points() const
Definition: fe.cc:1018
virtual FiniteElement< dim, spacedim > * clone() const =0
BlockIndices base_to_block_indices
Definition: fe.h:2242
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
virtual InternalDataBase * get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
bool has_support_points() const
Definition: fe.cc:972
types::global_dof_index first_block_of_base(const unsigned int b) const
Definition: fe.h:2852
static const unsigned int space_dimension
Definition: fe.h:582
bool has_face_support_points() const
Definition: fe.cc:1034
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
Definition: fe.cc:1236
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:2195
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2267
UpdateFlags update_each
Definition: fe.h:644
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2230
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:221
virtual InternalDataBase * get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature,::internal::FEValues::FiniteElementRelatedData< dim, spacedim > &output_data) const
Definition: fe.cc:1222