Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
fe.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe.h 31852 2013-12-03 13:54:24Z bangerth @f$
3 //
4 // Copyright (C) 1998 - 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_h
18 #define __deal2__fe_h
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/geometry_info.h>
22 #include <deal.II/fe/fe_base.h>
23 #include <deal.II/fe/fe_values_extractors.h>
24 #include <deal.II/fe/component_mask.h>
25 #include <deal.II/fe/block_mask.h>
26 
27 
29 
30 template <int dim, int spacedim> class FEValuesData;
31 template <int dim, int spacedim> class FEValuesBase;
32 template <int dim, int spacedim> class FEValues;
33 template <int dim, int spacedim> class FEFaceValues;
34 template <int dim, int spacedim> class FESubfaceValues;
35 template <int dim, int spacedim> class FESystem;
36 template <class POLY, int dim, int spacedim> class FE_PolyTensor;
37 namespace hp
38 {
39  template <int dim, int spacedim> class FECollection;
40 }
41 
42 
352 template <int dim, int spacedim=dim>
353 class FiniteElement : public Subscriptor,
354  public FiniteElementData<dim>
355 {
356 public:
366  class InternalDataBase : public Mapping<dim,spacedim>::InternalDataBase
367  {
368  public:
372  virtual ~InternalDataBase ();
373 
378  void initialize_2nd (const FiniteElement<dim,spacedim> *element,
379  const Mapping<dim,spacedim> &mapping,
380  const Quadrature<dim> &quadrature);
381 
389  std::vector<FEValues<dim,spacedim>*> differences;
390  };
391 
392 public:
396  FiniteElement (const FiniteElementData<dim> &fe_data,
397  const std::vector<bool> &restriction_is_additive_flags,
398  const std::vector<ComponentMask> &nonzero_components);
399 
404  virtual ~FiniteElement ();
405 
416  virtual std::string get_name () const = 0;
417 
439  const FiniteElement<dim,spacedim> &operator[] (const unsigned int fe_index) const;
440 
460  virtual double shape_value (const unsigned int i,
461  const Point<dim> &p) const;
462 
469  virtual double shape_value_component (const unsigned int i,
470  const Point<dim> &p,
471  const unsigned int component) const;
472 
489  virtual Tensor<1,dim> shape_grad (const unsigned int i,
490  const Point<dim> &p) const;
491 
498  virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
499  const Point<dim> &p,
500  const unsigned int component) const;
501 
517  virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
518  const Point<dim> &p) const;
519 
532  virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i,
533  const Point<dim> &p,
534  const unsigned int component) const;
545  virtual bool has_support_on_face (const unsigned int shape_index,
546  const unsigned int face_index) const;
547 
549 
573  virtual const FullMatrix<double> &
574  get_restriction_matrix (const unsigned int child,
576 
602  virtual const FullMatrix<double> &
603  get_prolongation_matrix (const unsigned int child,
605 
627  bool prolongation_is_implemented () const;
628 
645 
667  bool restriction_is_implemented () const;
668 
685 
686 
694  bool restriction_is_additive (const unsigned int index) const;
695 
710  const FullMatrix<double> &constraints (const ::internal::SubfaceCase<dim> &subface_case=::internal::SubfaceCase<dim>::case_isotropic) const;
711 
727  bool constraints_are_implemented (const ::internal::SubfaceCase<dim> &subface_case=::internal::SubfaceCase<dim>::case_isotropic) const;
728 
729 
751  virtual bool hp_constraints_are_implemented () const;
752 
753 
765  virtual void
767  FullMatrix<double> &matrix) const;
769 
787  virtual void
789  FullMatrix<double> &matrix) const;
790 
791 
803  virtual void
805  const unsigned int subface,
806  FullMatrix<double> &matrix) const;
808 
809 
832  virtual
833  std::vector<std::pair<unsigned int, unsigned int> >
835 
840  virtual
841  std::vector<std::pair<unsigned int, unsigned int> >
842  hp_line_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
843 
848  virtual
849  std::vector<std::pair<unsigned int, unsigned int> >
850  hp_quad_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
851 
860  virtual
863 
865 
875  bool operator == (const FiniteElement<dim,spacedim> &) const;
876 
907  std::pair<unsigned int, unsigned int>
908  system_to_component_index (const unsigned int index) const;
909 
919  unsigned int component_to_system_index(const unsigned int component,
920  const unsigned int index) const;
921 
931  std::pair<unsigned int, unsigned int>
932  face_system_to_component_index (const unsigned int index) const;
933 
986  virtual
987  unsigned int face_to_cell_index (const unsigned int face_dof_index,
988  const unsigned int face,
989  const bool face_orientation = true,
990  const bool face_flip = false,
991  const bool face_rotation = false) const;
992 
1001  unsigned int adjust_quad_dof_index_for_face_orientation (const unsigned int index,
1002  const bool face_orientation,
1003  const bool face_flip,
1004  const bool face_rotation) const;
1005 
1013  unsigned int adjust_line_dof_index_for_line_orientation (const unsigned int index,
1014  const bool line_orientation) const;
1015 
1032  const ComponentMask &
1033  get_nonzero_components (const unsigned int i) const;
1034 
1045  unsigned int
1046  n_nonzero_components (const unsigned int i) const;
1047 
1058  bool
1059  is_primitive (const unsigned int i) const;
1060 
1066 
1078  unsigned int n_base_elements () const;
1079 
1084  virtual
1086  base_element (const unsigned int index) const;
1087 
1094  unsigned int
1095  element_multiplicity (const unsigned int index) const;
1096 
1116  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1117  system_to_base_index (const unsigned int index) const;
1118 
1127  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1128  face_system_to_base_index (const unsigned int index) const;
1129 
1134  types::global_dof_index first_block_of_base (const unsigned int b) const;
1135 
1148  std::pair<unsigned int, unsigned int>
1149  component_to_base_index (const unsigned int component) const;
1150 
1151 
1156  std::pair<unsigned int,unsigned int>
1157  block_to_base_index (const unsigned int block) const;
1158 
1162  std::pair<unsigned int,types::global_dof_index>
1163  system_to_block_index (const unsigned int component) const;
1164 
1168  unsigned int
1169  component_to_block_index (const unsigned int component) const;
1170 
1172 
1190  component_mask (const FEValuesExtractors::Scalar &scalar) const;
1191 
1204  component_mask (const FEValuesExtractors::Vector &vector) const;
1205 
1218  component_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1219 
1233  component_mask (const BlockMask &block_mask) const;
1234 
1253  BlockMask
1254  block_mask (const FEValuesExtractors::Scalar &scalar) const;
1255 
1271  BlockMask
1272  block_mask (const FEValuesExtractors::Vector &vector) const;
1273 
1290  BlockMask
1291  block_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1292 
1313  BlockMask
1314  block_mask (const ComponentMask &component_mask) const;
1315 
1317 
1352  const std::vector<Point<dim> > &
1353  get_unit_support_points () const;
1354 
1369  bool has_support_points () const;
1370 
1383  virtual
1384  Point<dim>
1385  unit_support_point (const unsigned int index) const;
1386 
1414  const std::vector<Point<dim-1> > &
1416 
1425  bool has_face_support_points () const;
1426 
1431  virtual
1432  Point<dim-1>
1433  unit_face_support_point (const unsigned int index) const;
1434 
1441  const std::vector<Point<dim> > &
1443 
1451  bool has_generalized_support_points () const;
1452 
1456  const std::vector<Point<dim-1> > &
1457  get_generalized_face_support_points () const;
1458 
1468 
1478  virtual void interpolate(std::vector<double> &local_dofs,
1479  const std::vector<double> &values) const;
1480 
1490  virtual void interpolate(std::vector<double> &local_dofs,
1491  const std::vector<Vector<double> > &values,
1492  unsigned int offset = 0) const;
1493 
1498  virtual void interpolate(
1499  std::vector<double> &local_dofs,
1500  const VectorSlice<const std::vector<std::vector<double> > > &values) const;
1501 
1503 
1512  virtual std::size_t memory_consumption () const;
1518  DeclException1 (ExcShapeFunctionNotPrimitive,
1519  int,
1520  << "The shape function with index " << arg1
1521  << " is not primitive, i.e. it is vector-valued and "
1522  << "has more than one non-zero vector component. This "
1523  << "function cannot be called for these shape functions. "
1524  << "Maybe you want to use the same function with the "
1525  << "_component suffix?");
1531  DeclException0 (ExcFENotPrimitive);
1537  DeclException0 (ExcUnitShapeValuesDoNotExist);
1538 
1545  DeclException0 (ExcFEHasNoSupportPoints);
1546 
1553  DeclException0 (ExcEmbeddingVoid);
1554 
1562  DeclException0 (ExcProjectionVoid);
1563 
1571  DeclException0 (ExcConstraintsVoid);
1572 
1577  DeclException2 (ExcWrongInterfaceMatrixSize,
1578  int, int,
1579  << "The interface matrix has a size of " << arg1
1580  << "x" << arg2
1581  << ", which is not reasonable in the present dimension.");
1586  DeclException2 (ExcComponentIndexInvalid,
1587  int, int,
1588  << "The component-index pair (" << arg1 << ", " << arg2
1589  << ") is invalid, i.e. non-existent.");
1594  DeclException0 (ExcInterpolationNotImplemented);
1595 
1601  DeclException0 (ExcBoundaryFaceUsed);
1607  DeclException0 (ExcJacobiDeterminantHasWrongSign);
1608 
1609 protected:
1610 
1624  void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false,
1625  const bool isotropic_prolongation_only=false);
1626 
1639  std::vector<std::vector<FullMatrix<double> > > restriction;
1640 
1653  std::vector<std::vector<FullMatrix<double> > > prolongation;
1654 
1666 
1677  std::vector<Point<dim> > unit_support_points;
1678 
1684  std::vector<Point<dim-1> > unit_face_support_points;
1685 
1690  std::vector<Point<dim> > generalized_support_points;
1691 
1697 
1714 
1729 
1742  interface_constraints_size () const;
1743 
1747  void compute_2nd (const Mapping<dim,spacedim> &mapping,
1748  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
1749  const unsigned int offset,
1750  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
1751  InternalDataBase &fe_internal,
1752  FEValuesData<dim,spacedim> &data) const;
1753 
1759  static
1760  std::vector<unsigned int>
1761  compute_n_nonzero_components (const std::vector<ComponentMask> &nonzero_components);
1762 
1782  virtual UpdateFlags update_once (const UpdateFlags flags) const = 0;
1783 
1793  virtual UpdateFlags update_each (const UpdateFlags flags) const = 0;
1794 
1801  virtual FiniteElement<dim,spacedim> *clone() const = 0;
1802 
1803 private:
1807  std::vector< std::pair<unsigned int, unsigned int> > system_to_component_table;
1808 
1818  std::vector< std::pair<unsigned int, unsigned int> > face_system_to_component_table;
1819 
1836  std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
1838 
1842  std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
1844 
1850 
1871  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int> >
1873 
1909  const std::vector<bool> restriction_is_additive_flags;
1910 
1918  const std::vector<ComponentMask> nonzero_components;
1919 
1927  const std::vector<unsigned int> n_nonzero_components_table;
1928 
1934  static const double fd_step_length;
1935 
1942  virtual typename Mapping<dim,spacedim>::InternalDataBase *
1943  get_data (const UpdateFlags flags,
1944  const Mapping<dim,spacedim> &mapping,
1945  const Quadrature<dim> &quadrature) const = 0;
1946 
1953  virtual typename Mapping<dim,spacedim>::InternalDataBase *
1954  get_face_data (const UpdateFlags flags,
1955  const Mapping<dim,spacedim> &mapping,
1956  const Quadrature<dim-1> &quadrature) const;
1957 
1964  virtual typename Mapping<dim,spacedim>::InternalDataBase *
1965  get_subface_data (const UpdateFlags flags,
1966  const Mapping<dim,spacedim> &mapping,
1967  const Quadrature<dim-1> &quadrature) const;
1968 
1976  virtual void
1977  fill_fe_values (const Mapping<dim,spacedim> &mapping,
1978  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
1979  const Quadrature<dim> &quadrature,
1980  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
1981  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
1983  CellSimilarity::Similarity &cell_similarity) const = 0;
1984 
1992  virtual void
1994  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
1995  const unsigned int face_no,
1996  const Quadrature<dim-1> &quadrature,
1997  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
1998  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
1999  FEValuesData<dim,spacedim> &data) const = 0;
2000 
2008  virtual void
2010  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
2011  const unsigned int face_no,
2012  const unsigned int sub_no,
2013  const Quadrature<dim-1> &quadrature,
2014  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
2015  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
2016  FEValuesData<dim,spacedim> &data) const = 0;
2017 
2018  friend class InternalDataBase;
2019  friend class FEValuesBase<dim,spacedim>;
2020  friend class FEValues<dim,spacedim>;
2021  friend class FEFaceValues<dim,spacedim>;
2022  friend class FESubfaceValues<dim,spacedim>;
2023  template <int, int > friend class FESystem;
2024  template <class POLY, int dim_, int spacedim_> friend class FE_PolyTensor;
2025  friend class hp::FECollection<dim,spacedim>;
2026 
2027 };
2028 
2029 
2030 //----------------------------------------------------------------------//
2031 
2032 
2033 template <int dim, int spacedim>
2034 inline
2036 FiniteElement<dim,spacedim>::operator[] (const unsigned int fe_index) const
2037 {
2038  Assert (fe_index == 0,
2039  ExcMessage ("A fe_index of zero is the only index allowed here"));
2040  return *this;
2041 }
2042 
2043 
2044 
2045 template <int dim, int spacedim>
2046 inline
2047 std::pair<unsigned int,unsigned int>
2049 {
2050  Assert (index < system_to_component_table.size(),
2051  ExcIndexRange(index, 0, system_to_component_table.size()));
2052  Assert (is_primitive (index),
2054  return system_to_component_table[index];
2055 }
2056 
2057 
2058 
2059 template <int dim, int spacedim>
2060 inline
2061 unsigned int
2063 {
2064  return base_to_block_indices.size();
2065 }
2066 
2067 
2068 
2069 template <int dim, int spacedim>
2070 inline
2071 unsigned int
2073 {
2074  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
2075 }
2076 
2077 
2078 
2079 template <int dim, int spacedim>
2080 inline
2081 unsigned int
2083  const unsigned int index) const
2084 {
2085  std::vector< std::pair<unsigned int, unsigned int> >::const_iterator
2086  it = std::find(system_to_component_table.begin(), system_to_component_table.end(),
2087  std::pair<unsigned int, unsigned int>(component, index));
2088 
2089  Assert(it != system_to_component_table.end(), ExcComponentIndexInvalid(component, index));
2090  return std::distance(system_to_component_table.begin(), it);
2091 }
2092 
2093 
2094 
2095 template <int dim, int spacedim>
2096 inline
2097 std::pair<unsigned int,unsigned int>
2099 {
2100  Assert(index < face_system_to_component_table.size(),
2101  ExcIndexRange(index, 0, face_system_to_component_table.size()));
2102 
2103  // in debug mode, check whether the
2104  // function is primitive, since
2105  // otherwise the result may have no
2106  // meaning
2107  //
2108  // since the primitivity tables are
2109  // all geared towards cell dof
2110  // indices, rather than face dof
2111  // indices, we have to work a
2112  // little bit...
2113  //
2114  // in 1d, the face index is equal
2115  // to the cell index
2116  Assert (is_primitive(this->face_to_cell_index(index, 0)),
2118 
2119  return face_system_to_component_table[index];
2120 }
2121 
2122 
2123 
2124 template <int dim, int spacedim>
2125 inline
2126 std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
2128 {
2129  Assert (index < system_to_base_table.size(),
2130  ExcIndexRange(index, 0, system_to_base_table.size()));
2131  return system_to_base_table[index];
2132 }
2133 
2134 
2135 
2136 
2137 template <int dim, int spacedim>
2138 inline
2139 std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
2141 {
2142  Assert(index < face_system_to_base_table.size(),
2143  ExcIndexRange(index, 0, face_system_to_base_table.size()));
2144  return face_system_to_base_table[index];
2145 }
2146 
2147 
2148 
2149 template <int dim, int spacedim>
2150 inline
2153 {
2154  return base_to_block_indices.block_start(index);
2155 }
2156 
2157 
2158 
2159 template <int dim, int spacedim>
2160 inline
2161 std::pair<unsigned int,unsigned int>
2163 {
2164  Assert(index < component_to_base_table.size(),
2165  ExcIndexRange(index, 0, component_to_base_table.size()));
2166 
2167  return component_to_base_table[index].first;
2168 }
2169 
2170 
2171 
2172 template <int dim, int spacedim>
2173 inline
2174 std::pair<unsigned int,unsigned int>
2176 {
2177  return base_to_block_indices.global_to_local(index);
2178 }
2179 
2180 
2181 
2182 template <int dim, int spacedim>
2183 inline
2184 std::pair<unsigned int,types::global_dof_index>
2186 {
2187  Assert (index < this->dofs_per_cell,
2188  ExcIndexRange(index, 0, this->dofs_per_cell));
2189  // The block is computed simply as
2190  // first block of this base plus
2191  // the index within the base blocks
2192  return std::pair<unsigned int, types::global_dof_index>(
2193  first_block_of_base(system_to_base_table[index].first.first)
2194  + system_to_base_table[index].first.second,
2195  system_to_base_table[index].second);
2196 }
2197 
2198 
2199 
2200 template <int dim, int spacedim>
2201 inline
2202 bool
2204 {
2205  Assert(index < this->dofs_per_cell,
2206  ExcIndexRange(index, 0, this->dofs_per_cell));
2207  return restriction_is_additive_flags[index];
2208 }
2209 
2210 
2211 
2212 template <int dim, int spacedim>
2213 inline
2214 const ComponentMask &
2216 {
2217  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
2218  return nonzero_components[i];
2219 }
2220 
2221 
2222 
2223 template <int dim, int spacedim>
2224 inline
2225 unsigned int
2227 {
2228  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
2229  return n_nonzero_components_table[i];
2230 }
2231 
2232 
2233 
2234 template <int dim, int spacedim>
2235 inline
2236 bool
2237 FiniteElement<dim,spacedim>::is_primitive (const unsigned int i) const
2238 {
2239  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
2240 
2241  // return primitivity of a shape
2242  // function by checking whether it
2243  // has more than one non-zero
2244  // component or not. we could cache
2245  // this value in an array of bools,
2246  // but accessing a bit-vector (as
2247  // std::vector<bool> is) is
2248  // probably more expensive than
2249  // just comparing against 1
2250  //
2251  // for good measure, short circuit the test
2252  // if the entire FE is primitive
2253  return (is_primitive() ||
2254  (n_nonzero_components_table[i] == 1));
2255 }
2256 
2257 
2258 
2259 DEAL_II_NAMESPACE_CLOSE
2260 
2261 #endif
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index) const
Definition: fe.h:2098
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
virtual ~FiniteElement()
BlockIndices base_to_block_indices
Definition: fe.h:1849
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
std::vector< Point< dim > > unit_support_points
Definition: fe.h:1677
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
bool isotropic_restriction_is_implemented() const
virtual Point< dim-1 > unit_face_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_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const =0
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:1653
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_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
TableIndices< 2 > interface_constraints_size() const
Table< 2, int > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:1713
virtual void interpolate(std::vector< double > &local_dofs, const std::vector< double > &values) const
void initialize_2nd(const FiniteElement< dim, spacedim > *element, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature)
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:1909
std::vector< Point< dim-1 > > unit_face_support_points
Definition: fe.h:1684
::ExceptionBase & ExcMessage(std::string arg1)
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Auxiliary class aiding in the handling of block structures like in BlockVector or FESystem...
Definition: block_indices.h:55
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > face_system_to_base_table
Definition: fe.h:1843
const std::vector< Point< dim > > & get_unit_support_points() const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:1837
void compute_2nd(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int offset, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
unsigned int n_nonzero_components(const unsigned int i) const
Definition: fe.h:2226
bool is_primitive() const
std::vector< std::pair< unsigned int, unsigned int > > face_system_to_component_table
Definition: fe.h:1818
virtual FiniteElement< dim, spacedim > * clone() const =0
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:2048
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:1639
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 bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
Definition: fe.h:2185
virtual std::size_t memory_consumption() const
virtual Mapping< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature) const
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:1927
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
Definition: fe.h:2175
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > face_system_to_base_index(const unsigned int index) const
Definition: fe.h:2140
const std::vector< Point< dim-1 > > & get_unit_face_support_points() const
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:1728
Definition: fe.h:35
types::global_dof_index first_block_of_base(const unsigned int b) const
Definition: fe.h:2152
virtual bool hp_constraints_are_implemented() const
virtual std::string get_name() const =0
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
virtual UpdateFlags update_each(const UpdateFlags flags) const =0
const FiniteElement< dim, spacedim > & operator[](const unsigned int fe_index) const
Definition: fe.h:2036
unsigned int global_dof_index
Definition: types.h:100
bool prolongation_is_implemented() const
virtual UpdateFlags update_once(const UpdateFlags flags) const =0
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:299
DeclException0(ExcFENotPrimitive)
UpdateFlags
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
bool operator==(const FiniteElement< dim, spacedim > &) const
bool isotropic_prolongation_is_implemented() const
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:1918
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
Definition: fe.h:2127
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
unsigned int element_multiplicity(const unsigned int index) const
Definition: fe.h:2072
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:1690
virtual Mapping< dim, spacedim >::InternalDataBase * get_face_data(const UpdateFlags flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim-1 > &quadrature) const
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
Definition: fe.h:2082
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const
bool restriction_is_implemented() const
Definition: hp.h:103
virtual Point< dim > unit_support_point(const unsigned int index) const
bool has_support_points() const
bool has_generalized_face_support_points() const
std::vector< FEValues< dim, spacedim > * > differences
Definition: fe.h:389
bool restriction_is_additive(const unsigned int index) const
Definition: fe.h:2203
unsigned int component_to_block_index(const unsigned int component) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) 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_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data, CellSimilarity::Similarity &cell_similarity) const =0
bool has_generalized_support_points() const
unsigned int n_base_elements() const
Definition: fe.h:2062
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
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_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const =0
static const double fd_step_length
Definition: fe.h:1934
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const =0
DeclException2(ExcWrongInterfaceMatrixSize, int, int,<< "The interface matrix has a size of "<< arg1<< "x"<< arg2<< ", which is not reasonable in the present dimension.")
bool has_face_support_points() const
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:1872
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:1807
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
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
Definition: fe.h:2162
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
std::vector< Point< dim-1 > > generalized_face_support_points
Definition: fe.h:1696
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:2215
FullMatrix< double > interface_constraints
Definition: fe.h:1665
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
const std::vector< Point< dim > > & get_generalized_support_points() const