Reference documentation for deal.II version GIT 3a483e2de5 2022-11-28 08:40:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2022 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.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_fe_h
17 #define dealii_fe_h
18 
19 #include <deal.II/base/config.h>
20 
21 #include <deal.II/fe/block_mask.h>
23 #include <deal.II/fe/fe_data.h>
26 #include <deal.II/fe/mapping.h>
27 
29 #include <deal.II/lac/vector.h>
30 
31 #include <memory>
32 
33 
35 
36 // Forward declarations:
37 #ifndef DOXYGEN
38 template <int dim, int spacedim>
39 class FEValuesBase;
40 template <int dim, int spacedim>
41 class FEValues;
42 template <int dim, int spacedim>
43 class FEFaceValues;
44 template <int dim, int spacedim>
45 class FESubfaceValues;
46 namespace NonMatching
47 {
48  template <int dim>
49  class FEImmersedSurfaceValues;
50 }
51 template <int dim, int spacedim>
52 class FESystem;
53 #endif
54 
653 template <int dim, int spacedim = dim>
654 class FiniteElement : public Subscriptor, public FiniteElementData<dim>
655 {
656 public:
660  static constexpr unsigned int space_dimension = spacedim;
661 
687  {
688  public:
694 
698  virtual ~InternalDataBase() = default;
699 
704 
720 
724  virtual std::size_t
726  };
727 
728 public:
772  const std::vector<bool> & restriction_is_additive_flags,
773  const std::vector<ComponentMask> &nonzero_components);
774 
779 
784 
789  virtual ~FiniteElement() override = default;
790 
799  std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>
800  operator^(const unsigned int multiplicity) const;
801 
813  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
814  clone() const = 0;
815 
826  virtual std::string
827  get_name() const = 0;
828 
854  virtual double
855  shape_value(const unsigned int i, const Point<dim> &p) const;
856 
863  virtual double
864  shape_value_component(const unsigned int i,
865  const Point<dim> & p,
866  const unsigned int component) const;
867 
889  virtual Tensor<1, dim>
890  shape_grad(const unsigned int i, const Point<dim> &p) const;
891 
898  virtual Tensor<1, dim>
899  shape_grad_component(const unsigned int i,
900  const Point<dim> & p,
901  const unsigned int component) const;
902 
924  virtual Tensor<2, dim>
925  shape_grad_grad(const unsigned int i, const Point<dim> &p) const;
926 
933  virtual Tensor<2, dim>
934  shape_grad_grad_component(const unsigned int i,
935  const Point<dim> & p,
936  const unsigned int component) const;
937 
959  virtual Tensor<3, dim>
960  shape_3rd_derivative(const unsigned int i, const Point<dim> &p) const;
961 
968  virtual Tensor<3, dim>
969  shape_3rd_derivative_component(const unsigned int i,
970  const Point<dim> & p,
971  const unsigned int component) const;
972 
994  virtual Tensor<4, dim>
995  shape_4th_derivative(const unsigned int i, const Point<dim> &p) const;
996 
1003  virtual Tensor<4, dim>
1004  shape_4th_derivative_component(const unsigned int i,
1005  const Point<dim> & p,
1006  const unsigned int component) const;
1017  virtual bool
1018  has_support_on_face(const unsigned int shape_index,
1019  const unsigned int face_index) const;
1020 
1042  virtual const FullMatrix<double> &
1043  get_restriction_matrix(const unsigned int child,
1044  const RefinementCase<dim> &refinement_case =
1046 
1076  virtual const FullMatrix<double> &
1077  get_prolongation_matrix(const unsigned int child,
1078  const RefinementCase<dim> &refinement_case =
1080 
1102  bool
1104 
1120  bool
1122 
1144  bool
1146 
1162  bool
1164 
1165 
1174  bool
1175  restriction_is_additive(const unsigned int index) const;
1176 
1188  const FullMatrix<double> &
1189  constraints(const ::internal::SubfaceCase<dim> &subface_case =
1191 
1207  bool
1209  const ::internal::SubfaceCase<dim> &subface_case =
1211 
1212 
1234  virtual bool
1236 
1237 
1249  virtual void
1251  FullMatrix<double> & matrix) const;
1271  virtual void
1274  const unsigned int face_no = 0) const;
1275 
1276 
1288  virtual void
1290  const unsigned int subface,
1292  const unsigned int face_no = 0) const;
1316  virtual std::vector<std::pair<unsigned int, unsigned int>>
1318 
1323  virtual std::vector<std::pair<unsigned int, unsigned int>>
1325 
1330  virtual std::vector<std::pair<unsigned int, unsigned int>>
1332  const unsigned int face_no = 0) const;
1333 
1351  const unsigned int codim = 0) const;
1352 
1385  virtual bool
1387 
1392  bool
1394 
1430  std::pair<unsigned int, unsigned int>
1431  system_to_component_index(const unsigned int index) const;
1432 
1442  unsigned int
1443  component_to_system_index(const unsigned int component,
1444  const unsigned int index) const;
1445 
1455  std::pair<unsigned int, unsigned int>
1456  face_system_to_component_index(const unsigned int index,
1457  const unsigned int face_no = 0) const;
1458 
1467  unsigned int
1469  const unsigned int face_no,
1470  const bool face_orientation,
1471  const bool face_flip,
1472  const bool face_rotation) const;
1473 
1528  virtual unsigned int
1529  face_to_cell_index(const unsigned int face_dof_index,
1530  const unsigned int face,
1531  const bool face_orientation = true,
1532  const bool face_flip = false,
1533  const bool face_rotation = false) const;
1534 
1542  unsigned int
1544  const bool line_orientation) const;
1545 
1562  const ComponentMask &
1563  get_nonzero_components(const unsigned int i) const;
1564 
1575  unsigned int
1576  n_nonzero_components(const unsigned int i) const;
1577 
1586  bool
1587  is_primitive() const;
1588 
1599  bool
1600  is_primitive(const unsigned int i) const;
1601 
1613  unsigned int
1615 
1620  virtual const FiniteElement<dim, spacedim> &
1621  base_element(const unsigned int index) const;
1622 
1629  unsigned int
1630  element_multiplicity(const unsigned int index) const;
1631 
1705  get_sub_fe(const ComponentMask &mask) const;
1706 
1715  virtual const FiniteElement<dim, spacedim> &
1716  get_sub_fe(const unsigned int first_component,
1717  const unsigned int n_selected_components) const;
1718 
1741  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1742  system_to_base_index(const unsigned int index) const;
1743 
1752  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1753  face_system_to_base_index(const unsigned int index,
1754  const unsigned int face_no = 0) const;
1755 
1761  first_block_of_base(const unsigned int b) const;
1762 
1775  std::pair<unsigned int, unsigned int>
1776  component_to_base_index(const unsigned int component) const;
1777 
1778 
1783  std::pair<unsigned int, unsigned int>
1784  block_to_base_index(const unsigned int block) const;
1785 
1789  std::pair<unsigned int, types::global_dof_index>
1790  system_to_block_index(const unsigned int component) const;
1791 
1795  unsigned int
1796  component_to_block_index(const unsigned int component) const;
1797 
1819 
1834 
1850  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1851 
1868 
1889  BlockMask
1891 
1908  BlockMask
1910 
1928  BlockMask
1930 
1953  BlockMask
1955 
1971  virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
1973 
2011  const std::vector<Point<dim>> &
2013 
2031  bool
2033 
2046  virtual Point<dim>
2047  unit_support_point(const unsigned int index) const;
2048 
2075  const std::vector<Point<dim - 1>> &
2076  get_unit_face_support_points(const unsigned int face_no = 0) const;
2077 
2086  bool
2087  has_face_support_points(const unsigned int face_no = 0) const;
2088 
2093  virtual Point<dim - 1>
2094  unit_face_support_point(const unsigned int index,
2095  const unsigned int face_no = 0) const;
2096 
2124  const std::vector<Point<dim>> &
2126 
2136  bool
2138 
2183  get_associated_geometry_primitive(const unsigned int cell_dof_index) const;
2184 
2185 
2263  virtual void
2265  const std::vector<Vector<double>> &support_point_values,
2266  std::vector<double> & nodal_values) const;
2267 
2278  virtual std::size_t
2280 
2287  int,
2288  << "The shape function with index " << arg1
2289  << " is not primitive, i.e. it is vector-valued and "
2290  << "has more than one non-zero vector component. This "
2291  << "function cannot be called for these shape functions. "
2292  << "Maybe you want to use the same function with the "
2293  << "_component suffix?");
2307  "You are trying to access the values or derivatives of shape functions "
2308  "on the reference cell of an element that does not define its shape "
2309  "functions through mapping from the reference cell. Consequently, "
2310  "you cannot ask for shape function values or derivatives on the "
2311  "reference cell.");
2312 
2320  "You are trying to access the support points of a finite "
2321  "element that either has no support points at all, or for "
2322  "which the corresponding tables have not been implemented.");
2323 
2331  "You are trying to access the matrices that describe how "
2332  "to embed a finite element function on one cell into the "
2333  "finite element space on one of its children (i.e., the "
2334  "'embedding' or 'prolongation' matrices). However, the "
2335  "current finite element can either not define this sort of "
2336  "operation, or it has not yet been implemented.");
2337 
2346  "You are trying to access the matrices that describe how "
2347  "to restrict a finite element function from the children "
2348  "of one cell to the finite element space defined on their "
2349  "parent (i.e., the 'restriction' or 'projection' matrices). "
2350  "However, the current finite element can either not define "
2351  "this sort of operation, or it has not yet been "
2352  "implemented.");
2353 
2359  int,
2360  int,
2361  << "The interface matrix has a size of " << arg1 << 'x' << arg2
2362  << ", which is not reasonable for the current element "
2363  "in the present dimension.");
2369 
2370 protected:
2384  void
2386  const bool isotropic_restriction_only = false,
2387  const bool isotropic_prolongation_only = false);
2388 
2401  std::vector<std::vector<FullMatrix<double>>> restriction;
2402 
2415  std::vector<std::vector<FullMatrix<double>>> prolongation;
2416 
2428 
2439  std::vector<Point<dim>> unit_support_points;
2440 
2446  std::vector<std::vector<Point<dim - 1>>> unit_face_support_points;
2447 
2452  std::vector<Point<dim>> generalized_support_points;
2453 
2458  std::vector<std::vector<Point<dim - 1>>> generalized_face_support_points;
2459 
2476 
2491 
2495  std::vector<std::pair<unsigned int, unsigned int>> system_to_component_table;
2496 
2506  std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
2508 
2525  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2527 
2531  std::vector<
2532  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>>
2534 
2540 
2561  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2563 
2569  const std::vector<bool> restriction_is_additive_flags;
2570 
2578  const std::vector<ComponentMask> nonzero_components;
2579 
2587  const std::vector<unsigned int> n_nonzero_components_table;
2588 
2595 
2609 
2615  static std::vector<unsigned int>
2617  const std::vector<ComponentMask> &nonzero_components);
2618 
2639  virtual UpdateFlags
2640  requires_update_flags(const UpdateFlags update_flags) const = 0;
2641 
2718  virtual std::unique_ptr<InternalDataBase>
2719  get_data(const UpdateFlags update_flags,
2720  const Mapping<dim, spacedim> &mapping,
2721  const Quadrature<dim> & quadrature,
2723  FiniteElementRelatedData<dim, spacedim> &output_data) const = 0;
2724 
2766  virtual std::unique_ptr<InternalDataBase>
2767  get_face_data(const UpdateFlags update_flags,
2768  const Mapping<dim, spacedim> & mapping,
2769  const hp::QCollection<dim - 1> &quadrature,
2771  FiniteElementRelatedData<dim, spacedim> &output_data) const;
2772 
2776  virtual std::unique_ptr<InternalDataBase>
2778  const UpdateFlags update_flags,
2779  const Mapping<dim, spacedim> &mapping,
2780  const Quadrature<dim - 1> & quadrature,
2782  &output_data) const;
2783 
2825  virtual std::unique_ptr<InternalDataBase>
2827  const UpdateFlags update_flags,
2828  const Mapping<dim, spacedim> &mapping,
2829  const Quadrature<dim - 1> & quadrature,
2831  spacedim>
2832  &output_data) const;
2833 
2914  virtual void
2916  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2917  const CellSimilarity::Similarity cell_similarity,
2918  const Quadrature<dim> & quadrature,
2919  const Mapping<dim, spacedim> & mapping,
2920  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2921  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2922  spacedim>
2923  & mapping_data,
2924  const InternalDataBase &fe_internal,
2926  spacedim>
2927  &output_data) const = 0;
2928 
2971  virtual void
2973  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2974  const unsigned int face_no,
2975  const hp::QCollection<dim - 1> & quadrature,
2976  const Mapping<dim, spacedim> & mapping,
2977  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2978  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2979  spacedim>
2980  & mapping_data,
2981  const InternalDataBase &fe_internal,
2983  spacedim>
2984  &output_data) const;
2985 
2989  virtual void
2991  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2992  const unsigned int face_no,
2993  const Quadrature<dim - 1> & quadrature,
2994  const Mapping<dim, spacedim> & mapping,
2995  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2997  & mapping_data,
2998  const InternalDataBase &fe_internal,
3000  &output_data) const;
3001 
3047  virtual void
3049  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3050  const unsigned int face_no,
3051  const unsigned int sub_no,
3052  const Quadrature<dim - 1> & quadrature,
3053  const Mapping<dim, spacedim> & mapping,
3054  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3055  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
3056  spacedim>
3057  & mapping_data,
3058  const InternalDataBase &fe_internal,
3060  spacedim>
3061  &output_data) const = 0;
3062 
3063  friend class InternalDataBase;
3064  friend class FEValuesBase<dim, spacedim>;
3065  friend class FEValues<dim, spacedim>;
3066  friend class FEFaceValues<dim, spacedim>;
3067  friend class FESubfaceValues<dim, spacedim>;
3068  friend class NonMatching::FEImmersedSurfaceValues<dim>;
3069  friend class FESystem<dim, spacedim>;
3070 
3071  // explicitly check for sensible template arguments, but not on windows
3072  // because MSVC creates bogus warnings during normal compilation
3073 #ifndef DEAL_II_MSVC
3074  static_assert(dim <= spacedim,
3075  "The dimension <dim> of a FiniteElement must be less than or "
3076  "equal to the space dimension <spacedim> in which it lives.");
3077 #endif
3078 };
3079 
3080 
3081 //----------------------------------------------------------------------//
3082 
3083 
3084 template <int dim, int spacedim>
3085 inline std::pair<unsigned int, unsigned int>
3087  const unsigned int index) const
3088 {
3090  Assert(is_primitive(index),
3092  index)));
3093  return system_to_component_table[index];
3094 }
3095 
3096 
3097 
3098 template <int dim, int spacedim>
3099 inline unsigned int
3101 {
3102  return base_to_block_indices.size();
3103 }
3104 
3105 
3106 
3107 template <int dim, int spacedim>
3108 inline unsigned int
3110  const unsigned int index) const
3111 {
3112  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3113 }
3114 
3115 
3116 
3117 template <int dim, int spacedim>
3118 inline unsigned int
3120  const unsigned int component,
3121  const unsigned int index) const
3122 {
3123  AssertIndexRange(component, this->n_components());
3124  const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3125  std::find(system_to_component_table.begin(),
3127  std::pair<unsigned int, unsigned int>(component, index));
3128 
3129  Assert(it != system_to_component_table.end(),
3130  ExcMessage("You are asking for the number of the shape function "
3131  "within a system element that corresponds to vector "
3132  "component " +
3133  Utilities::int_to_string(component) +
3134  " and within this to "
3135  "index " +
3136  Utilities::int_to_string(index) +
3137  ". But no such "
3138  "shape function exists."));
3139  return std::distance(system_to_component_table.begin(), it);
3140 }
3141 
3142 
3143 
3144 template <int dim, int spacedim>
3145 inline std::pair<unsigned int, unsigned int>
3147  const unsigned int index,
3148  const unsigned int face_no) const
3149 {
3151  index,
3152  face_system_to_component_table[this->n_unique_faces() == 1 ? 0 : face_no]
3153  .size());
3154 
3155  // in debug mode, check whether the
3156  // function is primitive, since
3157  // otherwise the result may have no
3158  // meaning
3159  //
3160  // since the primitivity tables are
3161  // all geared towards cell dof
3162  // indices, rather than face dof
3163  // indices, we have to work a
3164  // little bit...
3165  //
3166  // in 1d, the face index is equal
3167  // to the cell index
3168  Assert(is_primitive(this->face_to_cell_index(index, face_no)),
3170  index)));
3171 
3172  return face_system_to_component_table[this->n_unique_faces() == 1 ?
3173  0 :
3174  face_no][index];
3175 }
3176 
3177 
3178 
3179 template <int dim, int spacedim>
3180 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3182  const unsigned int index) const
3183 {
3184  AssertIndexRange(index, system_to_base_table.size());
3185  return system_to_base_table[index];
3186 }
3187 
3188 
3189 
3190 template <int dim, int spacedim>
3191 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3193  const unsigned int index,
3194  const unsigned int face_no) const
3195 {
3197  index,
3198  face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3199  .size());
3200  return face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3201  [index];
3202 }
3203 
3204 
3205 
3206 template <int dim, int spacedim>
3209  const unsigned int index) const
3210 {
3211  return base_to_block_indices.block_start(index);
3212 }
3213 
3214 
3215 
3216 template <int dim, int spacedim>
3217 inline std::pair<unsigned int, unsigned int>
3219  const unsigned int index) const
3220 {
3222 
3223  return component_to_base_table[index].first;
3224 }
3225 
3226 
3227 
3228 template <int dim, int spacedim>
3229 inline std::pair<unsigned int, unsigned int>
3231  const unsigned int index) const
3232 {
3233  return base_to_block_indices.global_to_local(index);
3234 }
3235 
3236 
3237 
3238 template <int dim, int spacedim>
3239 inline std::pair<unsigned int, types::global_dof_index>
3241  const unsigned int index) const
3242 {
3243  AssertIndexRange(index, this->n_dofs_per_cell());
3244  // The block is computed simply as
3245  // first block of this base plus
3246  // the index within the base blocks
3247  return std::pair<unsigned int, types::global_dof_index>(
3249  system_to_base_table[index].first.second,
3250  system_to_base_table[index].second);
3251 }
3252 
3253 
3254 
3255 template <int dim, int spacedim>
3256 inline bool
3258  const unsigned int index) const
3259 {
3260  AssertIndexRange(index, this->n_dofs_per_cell());
3261  return restriction_is_additive_flags[index];
3262 }
3263 
3264 
3265 
3266 template <int dim, int spacedim>
3267 inline const ComponentMask &
3268 FiniteElement<dim, spacedim>::get_nonzero_components(const unsigned int i) const
3269 {
3270  AssertIndexRange(i, this->n_dofs_per_cell());
3271  return nonzero_components[i];
3272 }
3273 
3274 
3275 
3276 template <int dim, int spacedim>
3277 inline unsigned int
3278 FiniteElement<dim, spacedim>::n_nonzero_components(const unsigned int i) const
3279 {
3280  AssertIndexRange(i, this->n_dofs_per_cell());
3281  return n_nonzero_components_table[i];
3282 }
3283 
3284 
3285 
3286 template <int dim, int spacedim>
3287 inline bool
3289 {
3290  return cached_primitivity;
3291 }
3292 
3293 
3294 
3295 template <int dim, int spacedim>
3296 inline bool
3297 FiniteElement<dim, spacedim>::is_primitive(const unsigned int i) const
3298 {
3299  AssertIndexRange(i, this->n_dofs_per_cell());
3300 
3301  // return primitivity of a shape
3302  // function by checking whether it
3303  // has more than one non-zero
3304  // component or not. we could cache
3305  // this value in an array of bools,
3306  // but accessing a bit-vector (as
3307  // std::vector<bool> is) is
3308  // probably more expensive than
3309  // just comparing against 1
3310  //
3311  // for good measure, short circuit the test
3312  // if the entire FE is primitive
3313  return (is_primitive() || (n_nonzero_components_table[i] == 1));
3314 }
3315 
3316 
3317 
3318 template <int dim, int spacedim>
3319 inline GeometryPrimitive
3321  const unsigned int cell_dof_index) const
3322 {
3323  AssertIndexRange(cell_dof_index, this->n_dofs_per_cell());
3324 
3325  // just go through the usual cases, taking into account how DoFs
3326  // are enumerated on the reference cell
3327  if (cell_dof_index < this->get_first_line_index())
3329  else if (cell_dof_index < this->get_first_quad_index(0))
3330  return GeometryPrimitive::line;
3331  else if (cell_dof_index < this->get_first_hex_index())
3332  return GeometryPrimitive::quad;
3333  else
3334  return GeometryPrimitive::hex;
3335 }
3336 
3337 
3338 
3340 
3341 #endif
size_type block_size(const unsigned int i) const
unsigned int size() const
size_type block_start(const unsigned int i) const
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
unsigned int get_first_line_index() const
unsigned int n_dofs_per_cell() const
unsigned int get_first_quad_index(const unsigned int quad_no=0) const
unsigned int n_components() const
unsigned int n_unique_faces() const
unsigned int get_first_hex_index() const
virtual ~InternalDataBase()=default
InternalDataBase(const InternalDataBase &)=delete
UpdateFlags update_each
Definition: fe.h:719
virtual std::size_t memory_consumption() const
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
bool isotropic_prolongation_is_implemented() const
virtual std::string get_name() const =0
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2587
virtual std::unique_ptr< InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
const std::vector< Point< dim > > & get_generalized_support_points() const
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition: fe.h:2446
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
bool prolongation_is_implemented() const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
ComponentMask component_mask(const BlockMask &block_mask) const
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
ComponentMask component_mask(const FEValuesExtractors::SymmetricTensor< 2 > &sym_tensor) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const std::vector< Point< dim > > & get_unit_support_points() const
bool is_primitive(const unsigned int i) const
bool has_face_support_points(const unsigned int face_no=0) const
const bool cached_primitivity
Definition: fe.h:2594
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
bool operator!=(const FiniteElement< dim, spacedim > &) const
virtual std::unique_ptr< InternalDataBase > get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
ComponentMask component_mask(const FEValuesExtractors::Vector &vector) const
bool has_support_points() const
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2401
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
virtual const FiniteElement< dim, spacedim > & get_sub_fe(const unsigned int first_component, const unsigned int n_selected_components) const
TableIndices< 2 > interface_constraints_size() const
bool has_generalized_support_points() const
std::vector< Table< 2, int > > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2475
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double >> &support_point_values, std::vector< double > &nodal_values) const
const ComponentMask & get_nonzero_components(const unsigned int i) const
static constexpr unsigned int space_dimension
Definition: fe.h:660
BlockIndices base_to_block_indices
Definition: fe.h:2539
bool is_primitive() const
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const =0
bool isotropic_restriction_is_implemented() const
std::vector< std::vector< Point< dim - 1 > > > generalized_face_support_points
Definition: fe.h:2458
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:2490
virtual Point< dim > unit_support_point(const unsigned int index) const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > face_system_to_base_index(const unsigned int index, const unsigned int face_no=0) const
bool restriction_is_additive(const unsigned int index) const
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
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > operator^(const unsigned int multiplicity) const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
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::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
FiniteElement(const FiniteElement< dim, spacedim > &)=default
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const
unsigned int component_to_block_index(const unsigned int component) const
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
virtual std::unique_ptr< InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const =0
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
virtual std::unique_ptr< InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2526
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
BlockMask block_mask(const FEValuesExtractors::Vector &vector) const
BlockMask block_mask(const ComponentMask &component_mask) const
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2495
const FiniteElement< dim, spacedim > & get_sub_fe(const ComponentMask &mask) const
unsigned int element_multiplicity(const unsigned int index) const
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2569
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index, const unsigned int face_no=0) const
virtual ~FiniteElement() override=default
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2562
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > face_system_to_base_table
Definition: fe.h:2533
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2439
bool restriction_is_implemented() const
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > face_system_to_component_table
Definition: fe.h:2507
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::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
FiniteElement(FiniteElement< dim, spacedim > &&)=default
unsigned int n_nonzero_components(const unsigned int i) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) const
const std::vector< Point< dim - 1 > > & get_unit_face_support_points(const unsigned int face_no=0) const
FullMatrix< double > interface_constraints
Definition: fe.h:2427
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index, const unsigned int face_no=0) const
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
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::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
BlockMask block_mask(const FEValuesExtractors::SymmetricTensor< 2 > &sym_tensor) const
unsigned int n_base_elements() const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const
virtual std::size_t memory_consumption() const
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2452
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
unsigned int adjust_quad_dof_index_for_face_orientation(const unsigned int index, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation) const
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2578
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual bool hp_constraints_are_implemented() const
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2415
types::global_dof_index first_block_of_base(const unsigned int b) const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
Definition: point.h:111
Definition: vector.h:109
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
Point< 2 > first
Definition: grid_out.cc:4605
static ::ExceptionBase & ExcFEHasNoSupportPoints()
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
static ::ExceptionBase & ExcEmbeddingVoid()
#define Assert(cond, exc)
Definition: exceptions.h:1501
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
#define AssertIndexRange(index, range)
Definition: exceptions.h:1760
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcInterpolationNotImplemented()
static ::ExceptionBase & ExcProjectionVoid()
static ::ExceptionBase & ExcMessage(std::string arg1)
UpdateFlags
@ matrix
Contents is actually a matrix.
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:473
unsigned int global_dof_index
Definition: types.h:81