Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 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>
22 #include <deal.II/fe/component_mask.h>
23 #include <deal.II/fe/fe_base.h>
24 #include <deal.II/fe/fe_update_flags.h>
25 #include <deal.II/fe/fe_values_extractors.h>
26 #include <deal.II/fe/mapping.h>
27 
28 #include <deal.II/lac/full_matrix.h>
29 
30 #include <memory>
31 
32 
33 DEAL_II_NAMESPACE_OPEN
34 
35 template <int dim, int spacedim>
37 template <int dim, int spacedim>
38 class FEValues;
39 template <int dim, int spacedim>
41 template <int dim, int spacedim>
43 template <int dim, int spacedim>
44 class FESystem;
45 
647 template <int dim, int spacedim = dim>
648 class FiniteElement : public Subscriptor, public FiniteElementData<dim>
649 {
650 public:
654  static const unsigned int space_dimension = spacedim;
655 
683  {
684  public:
690 
694  virtual ~InternalDataBase() = default;
695 
699  InternalDataBase(const InternalDataBase &) = delete;
700 
716 
720  virtual std::size_t
721  memory_consumption() const;
722  };
723 
724 public:
767  FiniteElement(const FiniteElementData<dim> & fe_data,
768  const std::vector<bool> & restriction_is_additive_flags,
769  const std::vector<ComponentMask> &nonzero_components);
770 
774  FiniteElement(FiniteElement<dim, spacedim> &&) = default; // NOLINT
775 
779  FiniteElement(const FiniteElement<dim, spacedim> &) = default;
780 
785  virtual ~FiniteElement() override = default;
786 
795  std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>
796  operator^(const unsigned int multiplicity) const;
797 
809  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
810  clone() const = 0;
811 
822  virtual std::string
823  get_name() const = 0;
824 
847  operator[](const unsigned int fe_index) const;
848 
874  virtual double
875  shape_value(const unsigned int i, const Point<dim> &p) const;
876 
883  virtual double
884  shape_value_component(const unsigned int i,
885  const Point<dim> & p,
886  const unsigned int component) const;
887 
909  virtual Tensor<1, dim>
910  shape_grad(const unsigned int i, const Point<dim> &p) const;
911 
918  virtual Tensor<1, dim>
919  shape_grad_component(const unsigned int i,
920  const Point<dim> & p,
921  const unsigned int component) const;
922 
944  virtual Tensor<2, dim>
945  shape_grad_grad(const unsigned int i, const Point<dim> &p) const;
946 
953  virtual Tensor<2, dim>
954  shape_grad_grad_component(const unsigned int i,
955  const Point<dim> & p,
956  const unsigned int component) const;
957 
979  virtual Tensor<3, dim>
980  shape_3rd_derivative(const unsigned int i, const Point<dim> &p) const;
981 
988  virtual Tensor<3, dim>
989  shape_3rd_derivative_component(const unsigned int i,
990  const Point<dim> & p,
991  const unsigned int component) const;
992 
1014  virtual Tensor<4, dim>
1015  shape_4th_derivative(const unsigned int i, const Point<dim> &p) const;
1016 
1023  virtual Tensor<4, dim>
1024  shape_4th_derivative_component(const unsigned int i,
1025  const Point<dim> & p,
1026  const unsigned int component) const;
1037  virtual bool
1038  has_support_on_face(const unsigned int shape_index,
1039  const unsigned int face_index) const;
1040 
1042 
1062  virtual const FullMatrix<double> &
1063  get_restriction_matrix(const unsigned int child,
1064  const RefinementCase<dim> &refinement_case =
1066 
1096  virtual const FullMatrix<double> &
1097  get_prolongation_matrix(const unsigned int child,
1098  const RefinementCase<dim> &refinement_case =
1100 
1122  bool
1124 
1140  bool
1142 
1164  bool
1166 
1182  bool
1184 
1185 
1194  bool
1195  restriction_is_additive(const unsigned int index) const;
1196 
1208  const FullMatrix<double> &
1209  constraints(const ::internal::SubfaceCase<dim> &subface_case =
1211 
1227  bool
1229  const ::internal::SubfaceCase<dim> &subface_case =
1231 
1232 
1254  virtual bool
1256 
1257 
1269  virtual void
1271  FullMatrix<double> & matrix) const;
1273 
1291  virtual void
1293  FullMatrix<double> &matrix) const;
1294 
1295 
1307  virtual void
1309  const unsigned int subface,
1310  FullMatrix<double> &matrix) const;
1312 
1313 
1334  virtual std::vector<std::pair<unsigned int, unsigned int>>
1336 
1341  virtual std::vector<std::pair<unsigned int, unsigned int>>
1343 
1348  virtual std::vector<std::pair<unsigned int, unsigned int>>
1350 
1362  DEAL_II_DEPRECATED virtual FiniteElementDomination::Domination
1364  const FiniteElement<dim, spacedim> &fe_other) const final;
1365 
1383  const unsigned int codim = 0) const;
1384 
1386 
1417  virtual bool
1418  operator==(const FiniteElement<dim, spacedim> &fe) const;
1419 
1424  bool
1426 
1462  std::pair<unsigned int, unsigned int>
1463  system_to_component_index(const unsigned int index) const;
1464 
1474  unsigned int
1475  component_to_system_index(const unsigned int component,
1476  const unsigned int index) const;
1477 
1487  std::pair<unsigned int, unsigned int>
1488  face_system_to_component_index(const unsigned int index) const;
1489 
1498  unsigned int
1499  adjust_quad_dof_index_for_face_orientation(const unsigned int index,
1500  const bool face_orientation,
1501  const bool face_flip,
1502  const bool face_rotation) const;
1503 
1558  virtual unsigned int
1559  face_to_cell_index(const unsigned int face_dof_index,
1560  const unsigned int face,
1561  const bool face_orientation = true,
1562  const bool face_flip = false,
1563  const bool face_rotation = false) const;
1564 
1572  unsigned int
1573  adjust_line_dof_index_for_line_orientation(const unsigned int index,
1574  const bool line_orientation) const;
1575 
1592  const ComponentMask &
1593  get_nonzero_components(const unsigned int i) const;
1594 
1605  unsigned int
1606  n_nonzero_components(const unsigned int i) const;
1607 
1616  bool
1617  is_primitive() const;
1618 
1629  bool
1630  is_primitive(const unsigned int i) const;
1631 
1643  unsigned int
1644  n_base_elements() const;
1645 
1650  virtual const FiniteElement<dim, spacedim> &
1651  base_element(const unsigned int index) const;
1652 
1659  unsigned int
1660  element_multiplicity(const unsigned int index) const;
1661 
1735  get_sub_fe(const ComponentMask &mask) const;
1736 
1745  virtual const FiniteElement<dim, spacedim> &
1746  get_sub_fe(const unsigned int first_component,
1747  const unsigned int n_selected_components) const;
1748 
1771  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1772  system_to_base_index(const unsigned int index) const;
1773 
1782  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1783  face_system_to_base_index(const unsigned int index) const;
1784 
1790  first_block_of_base(const unsigned int b) const;
1791 
1804  std::pair<unsigned int, unsigned int>
1805  component_to_base_index(const unsigned int component) const;
1806 
1807 
1812  std::pair<unsigned int, unsigned int>
1813  block_to_base_index(const unsigned int block) const;
1814 
1818  std::pair<unsigned int, types::global_dof_index>
1819  system_to_block_index(const unsigned int component) const;
1820 
1824  unsigned int
1825  component_to_block_index(const unsigned int component) const;
1826 
1828 
1847  component_mask(const FEValuesExtractors::Scalar &scalar) const;
1848 
1862  component_mask(const FEValuesExtractors::Vector &vector) const;
1863 
1879  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1880 
1896  component_mask(const BlockMask &block_mask) const;
1897 
1918  BlockMask
1919  block_mask(const FEValuesExtractors::Scalar &scalar) const;
1920 
1937  BlockMask
1938  block_mask(const FEValuesExtractors::Vector &vector) const;
1939 
1957  BlockMask
1958  block_mask(const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1959 
1982  BlockMask
1983  block_mask(const ComponentMask &component_mask) const;
1984 
2000  virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
2001  get_constant_modes() const;
2002 
2004 
2040  const std::vector<Point<dim>> &
2041  get_unit_support_points() const;
2042 
2060  bool
2061  has_support_points() const;
2062 
2075  virtual Point<dim>
2076  unit_support_point(const unsigned int index) const;
2077 
2104  const std::vector<Point<dim - 1>> &
2106 
2115  bool
2116  has_face_support_points() const;
2117 
2122  virtual Point<dim - 1>
2123  unit_face_support_point(const unsigned int index) const;
2124 
2137  const std::vector<Point<dim>> &
2139 
2149  bool
2151 
2160  DEAL_II_DEPRECATED
2161  const std::vector<Point<dim - 1>> &
2163 
2176  DEAL_II_DEPRECATED
2177  bool
2179 
2224  get_associated_geometry_primitive(const unsigned int cell_dof_index) const;
2225 
2226 
2304  virtual void
2306  const std::vector<Vector<double>> &support_point_values,
2307  std::vector<double> & nodal_values) const;
2308 
2310 
2319  virtual std::size_t
2320  memory_consumption() const;
2321 
2328  int,
2329  << "The shape function with index " << arg1
2330  << " is not primitive, i.e. it is vector-valued and "
2331  << "has more than one non-zero vector component. This "
2332  << "function cannot be called for these shape functions. "
2333  << "Maybe you want to use the same function with the "
2334  << "_component suffix?");
2348  "You are trying to access the values or derivatives of shape functions "
2349  "on the reference cell of an element that does not define its shape "
2350  "functions through mapping from the reference cell. Consequently, "
2351  "you cannot ask for shape function values or derivatives on the "
2352  "reference cell.");
2353 
2361  "You are trying to access the support points of a finite "
2362  "element that either has no support points at all, or for "
2363  "which the corresponding tables have not been implemented.");
2364 
2372  "You are trying to access the matrices that describe how "
2373  "to embed a finite element function on one cell into the "
2374  "finite element space on one of its children (i.e., the "
2375  "'embedding' or 'prolongation' matrices). However, the "
2376  "current finite element can either not define this sort of "
2377  "operation, or it has not yet been implemented.");
2378 
2387  "You are trying to access the matrices that describe how "
2388  "to restrict a finite element function from the children "
2389  "of one cell to the finite element space defined on their "
2390  "parent (i.e., the 'restriction' or 'projection' matrices). "
2391  "However, the current finite element can either not define "
2392  "this sort of operation, or it has not yet been "
2393  "implemented.");
2394 
2400  int,
2401  int,
2402  << "The interface matrix has a size of " << arg1 << "x" << arg2
2403  << ", which is not reasonable for the current element "
2404  "in the present dimension.");
2410 
2411 protected:
2425  void
2427  const bool isotropic_restriction_only = false,
2428  const bool isotropic_prolongation_only = false);
2429 
2442  std::vector<std::vector<FullMatrix<double>>> restriction;
2443 
2456  std::vector<std::vector<FullMatrix<double>>> prolongation;
2457 
2469 
2480  std::vector<Point<dim>> unit_support_points;
2481 
2487  std::vector<Point<dim - 1>> unit_face_support_points;
2488 
2493  std::vector<Point<dim>> generalized_support_points;
2494 
2500 
2517 
2532 
2536  std::vector<std::pair<unsigned int, unsigned int>> system_to_component_table;
2537 
2547  std::vector<std::pair<unsigned int, unsigned int>>
2549 
2566  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2568 
2572  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2574 
2580 
2601  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2603 
2609  const std::vector<bool> restriction_is_additive_flags;
2610 
2618  const std::vector<ComponentMask> nonzero_components;
2619 
2627  const std::vector<unsigned int> n_nonzero_components_table;
2628 
2635 
2649 
2655  static std::vector<unsigned int>
2657  const std::vector<ComponentMask> &nonzero_components);
2658 
2679  virtual UpdateFlags
2680  requires_update_flags(const UpdateFlags update_flags) const = 0;
2681 
2758  virtual std::unique_ptr<InternalDataBase>
2759  get_data(const UpdateFlags update_flags,
2760  const Mapping<dim, spacedim> &mapping,
2761  const Quadrature<dim> & quadrature,
2763  FiniteElementRelatedData<dim, spacedim> &output_data) const = 0;
2764 
2806  virtual std::unique_ptr<InternalDataBase>
2807  get_face_data(const UpdateFlags update_flags,
2808  const Mapping<dim, spacedim> &mapping,
2809  const Quadrature<dim - 1> & quadrature,
2811  FiniteElementRelatedData<dim, spacedim> &output_data) const;
2812 
2854  virtual std::unique_ptr<InternalDataBase>
2856  const UpdateFlags update_flags,
2857  const Mapping<dim, spacedim> &mapping,
2858  const Quadrature<dim - 1> & quadrature,
2860  spacedim>
2861  &output_data) const;
2862 
2943  virtual void
2945  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2946  const CellSimilarity::Similarity cell_similarity,
2947  const Quadrature<dim> & quadrature,
2948  const Mapping<dim, spacedim> & mapping,
2949  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2950  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2951  spacedim>
2952  & mapping_data,
2953  const InternalDataBase &fe_internal,
2955  spacedim>
2956  &output_data) const = 0;
2957 
3000  virtual void
3002  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3003  const unsigned int face_no,
3004  const Quadrature<dim - 1> & quadrature,
3005  const Mapping<dim, spacedim> & mapping,
3006  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3007  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
3008  spacedim>
3009  & mapping_data,
3010  const InternalDataBase &fe_internal,
3012  spacedim>
3013  &output_data) const = 0;
3014 
3060  virtual void
3062  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3063  const unsigned int face_no,
3064  const unsigned int sub_no,
3065  const Quadrature<dim - 1> & quadrature,
3066  const Mapping<dim, spacedim> & mapping,
3067  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3068  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
3069  spacedim>
3070  & mapping_data,
3071  const InternalDataBase &fe_internal,
3073  spacedim>
3074  &output_data) const = 0;
3075 
3076  friend class InternalDataBase;
3077  friend class FEValuesBase<dim, spacedim>;
3078  friend class FEValues<dim, spacedim>;
3079  friend class FEFaceValues<dim, spacedim>;
3080  friend class FESubfaceValues<dim, spacedim>;
3081  friend class FESystem<dim, spacedim>;
3082 
3083  // explicitly check for sensible template arguments, but not on windows
3084  // because MSVC creates bogus warnings during normal compilation
3085 #ifndef DEAL_II_MSVC
3086  static_assert(dim <= spacedim,
3087  "The dimension <dim> of a FiniteElement must be less than or "
3088  "equal to the space dimension <spacedim> in which it lives.");
3089 #endif
3090 };
3091 
3092 
3093 //----------------------------------------------------------------------//
3094 
3095 
3096 template <int dim, int spacedim>
3098  operator[](const unsigned int fe_index) const
3099 {
3100  (void)fe_index;
3101  Assert(fe_index == 0,
3102  ExcMessage("A fe_index of zero is the only index allowed here"));
3103  return *this;
3104 }
3105 
3106 
3107 
3108 template <int dim, int spacedim>
3109 inline std::pair<unsigned int, unsigned int>
3111  const unsigned int index) const
3112 {
3113  Assert(index < system_to_component_table.size(),
3114  ExcIndexRange(index, 0, system_to_component_table.size()));
3115  Assert(is_primitive(index),
3117  index)));
3118  return system_to_component_table[index];
3119 }
3120 
3121 
3122 
3123 template <int dim, int spacedim>
3124 inline unsigned int
3126 {
3127  return base_to_block_indices.size();
3128 }
3129 
3130 
3131 
3132 template <int dim, int spacedim>
3133 inline unsigned int
3135  const unsigned int index) const
3136 {
3137  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3138 }
3139 
3140 
3141 
3142 template <int dim, int spacedim>
3143 inline unsigned int
3145  const unsigned int component,
3146  const unsigned int index) const
3147 {
3148  AssertIndexRange(component, this->n_components());
3149  const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3150  std::find(system_to_component_table.begin(),
3152  std::pair<unsigned int, unsigned int>(component, index));
3153 
3154  Assert(it != system_to_component_table.end(),
3155  ExcMessage("You are asking for the number of the shape function "
3156  "within a system element that corresponds to vector "
3157  "component " +
3158  Utilities::int_to_string(component) +
3159  " and within this to "
3160  "index " +
3161  Utilities::int_to_string(index) +
3162  ". But no such "
3163  "shape function exists."));
3164  return std::distance(system_to_component_table.begin(), it);
3165 }
3166 
3167 
3168 
3169 template <int dim, int spacedim>
3170 inline std::pair<unsigned int, unsigned int>
3172  const unsigned int index) const
3173 {
3174  Assert(index < face_system_to_component_table.size(),
3175  ExcIndexRange(index, 0, face_system_to_component_table.size()));
3176 
3177  // in debug mode, check whether the
3178  // function is primitive, since
3179  // otherwise the result may have no
3180  // meaning
3181  //
3182  // since the primitivity tables are
3183  // all geared towards cell dof
3184  // indices, rather than face dof
3185  // indices, we have to work a
3186  // little bit...
3187  //
3188  // in 1d, the face index is equal
3189  // to the cell index
3190  Assert(is_primitive(this->face_to_cell_index(index, 0)),
3192  index)));
3193 
3194  return face_system_to_component_table[index];
3195 }
3196 
3197 
3198 
3199 template <int dim, int spacedim>
3200 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3202  const unsigned int index) const
3203 {
3204  Assert(index < system_to_base_table.size(),
3205  ExcIndexRange(index, 0, system_to_base_table.size()));
3206  return system_to_base_table[index];
3207 }
3208 
3209 
3210 
3211 template <int dim, int spacedim>
3212 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3214  const unsigned int index) const
3215 {
3216  Assert(index < face_system_to_base_table.size(),
3217  ExcIndexRange(index, 0, face_system_to_base_table.size()));
3218  return face_system_to_base_table[index];
3219 }
3220 
3221 
3222 
3223 template <int dim, int spacedim>
3226  const unsigned int index) const
3227 {
3228  return base_to_block_indices.block_start(index);
3229 }
3230 
3231 
3232 
3233 template <int dim, int spacedim>
3234 inline std::pair<unsigned int, unsigned int>
3236  const unsigned int index) const
3237 {
3238  Assert(index < component_to_base_table.size(),
3239  ExcIndexRange(index, 0, component_to_base_table.size()));
3240 
3241  return component_to_base_table[index].first;
3242 }
3243 
3244 
3245 
3246 template <int dim, int spacedim>
3247 inline std::pair<unsigned int, unsigned int>
3249  const unsigned int index) const
3250 {
3251  return base_to_block_indices.global_to_local(index);
3252 }
3253 
3254 
3255 
3256 template <int dim, int spacedim>
3257 inline std::pair<unsigned int, types::global_dof_index>
3259  const unsigned int index) const
3260 {
3261  Assert(index < this->dofs_per_cell,
3262  ExcIndexRange(index, 0, this->dofs_per_cell));
3263  // The block is computed simply as
3264  // first block of this base plus
3265  // the index within the base blocks
3266  return std::pair<unsigned int, types::global_dof_index>(
3267  first_block_of_base(system_to_base_table[index].first.first) +
3268  system_to_base_table[index].first.second,
3269  system_to_base_table[index].second);
3270 }
3271 
3272 
3273 
3274 template <int dim, int spacedim>
3275 inline bool
3277  const unsigned int index) const
3278 {
3279  Assert(index < this->dofs_per_cell,
3280  ExcIndexRange(index, 0, this->dofs_per_cell));
3281  return restriction_is_additive_flags[index];
3282 }
3283 
3284 
3285 
3286 template <int dim, int spacedim>
3287 inline const ComponentMask &
3289 {
3290  Assert(i < this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
3291  return nonzero_components[i];
3292 }
3293 
3294 
3295 
3296 template <int dim, int spacedim>
3297 inline unsigned int
3299 {
3300  Assert(i < this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
3301  return n_nonzero_components_table[i];
3302 }
3303 
3304 
3305 
3306 template <int dim, int spacedim>
3307 inline bool
3309 {
3310  return cached_primitivity;
3311 }
3312 
3313 
3314 
3315 template <int dim, int spacedim>
3316 inline bool
3318 {
3319  Assert(i < this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
3320 
3321  // return primitivity of a shape
3322  // function by checking whether it
3323  // has more than one non-zero
3324  // component or not. we could cache
3325  // this value in an array of bools,
3326  // but accessing a bit-vector (as
3327  // std::vector<bool> is) is
3328  // probably more expensive than
3329  // just comparing against 1
3330  //
3331  // for good measure, short circuit the test
3332  // if the entire FE is primitive
3333  return (is_primitive() || (n_nonzero_components_table[i] == 1));
3334 }
3335 
3336 
3337 
3338 template <int dim, int spacedim>
3339 inline GeometryPrimitive
3341  const unsigned int cell_dof_index) const
3342 {
3343  Assert(cell_dof_index < this->dofs_per_cell,
3344  ExcIndexRange(cell_dof_index, 0, this->dofs_per_cell));
3345 
3346  // just go through the usual cases, taking into account how DoFs
3347  // are enumerated on the reference cell
3348  if (cell_dof_index < this->first_line_index)
3350  else if (cell_dof_index < this->first_quad_index)
3351  return GeometryPrimitive::line;
3352  else if (cell_dof_index < this->first_hex_index)
3353  return GeometryPrimitive::quad;
3354  else
3355  return GeometryPrimitive::hex;
3356 }
3357 
3358 
3359 
3360 DEAL_II_NAMESPACE_CLOSE
3361 
3362 #endif
static ::ExceptionBase & ExcFEHasNoSupportPoints()
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
const unsigned int first_hex_index
Definition: fe_base.h:258
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const =0
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:465
bool restriction_is_additive(const unsigned int index) const
Definition: fe.h:3276
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:216
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2442
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:541
virtual std::size_t memory_consumption() const
Definition: fe.cc:49
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:896
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
Definition: fe.cc:985
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index) const
Definition: fe.cc:1118
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2493
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:229
FullMatrix< double > interface_constraints
Definition: fe.h:2468
bool operator!=(const FiniteElement< dim, spacedim > &) const
Definition: fe.cc:1000
unsigned int n_nonzero_components(const unsigned int i) const
Definition: fe.h:3298
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:57
const FiniteElement< dim, spacedim > & get_sub_fe(const ComponentMask &mask) const
Definition: fe.cc:1142
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:170
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:880
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
virtual Point< dim > unit_support_point(const unsigned int index) const
Definition: fe.cc:1056
const std::vector< Point< dim - 1 > > & get_unit_face_support_points() const
Definition: fe.cc:1069
bool restriction_is_implemented() const
Definition: fe.cc:724
const std::vector< Point< dim > > & get_generalized_support_points() const
Definition: fe.cc:1034
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:182
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:240
#define AssertIndexRange(index, range)
Definition: exceptions.h:1637
virtual ~InternalDataBase()=default
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:951
bool isotropic_prolongation_is_implemented() const
Definition: fe.cc:752
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2602
bool has_generalized_support_points() const
Definition: fe.cc:1047
std::vector< Point< dim - 1 > > unit_face_support_points
Definition: fe.h:2487
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:379
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const =0
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:1245
size_type block_size(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:335
bool is_primitive() const
Definition: fe.h:3308
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 override
Definition: fe_system.cc:869
static ::ExceptionBase & ExcInterpolationNotImplemented()
std::vector< Point< dim - 1 > > generalized_face_support_points
Definition: fe.h:2499
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2627
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
bool has_face_support_points() const
Definition: fe.cc:1085
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2480
bool has_support_points() const
Definition: fe.cc:1025
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
Definition: fe.cc:1278
unsigned int component_to_block_index(const unsigned int component) const
Definition: fe.cc:366
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcMessage(std::string arg1)
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:830
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:193
const unsigned int first_quad_index
Definition: fe_base.h:253
virtual std::size_t memory_consumption() const
Definition: fe.cc:1224
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2456
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
Definition: fe.h:3235
#define Assert(cond, exc)
Definition: exceptions.h:1407
unsigned int element_multiplicity(const unsigned int index) const
Definition: fe.h:3134
UpdateFlags
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
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const final
Definition: fe.cc:962
Abstract base class for mapping classes.
Definition: dof_tools.h:57
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
Definition: fe.h:3340
Definition: fe.h:44
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:3288
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
#define DeclException0(Exception0)
Definition: exceptions.h:473
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
Definition: fe.cc:912
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:808
static ::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index) const
Definition: fe.h:3171
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 =0
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
Definition: fe.cc:1296
virtual std::string get_name() const =0
static ::ExceptionBase & ExcProjectionVoid()
bool has_generalized_face_support_points() const
Definition: fe.cc:1109
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:929
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:551
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:383
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
Definition: fe.h:3248
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double >> &support_point_values, std::vector< double > &nodal_values) const
Definition: fe.cc:1207
const unsigned int dofs_per_cell
Definition: fe_base.h:282
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
types::global_dof_index first_block_of_base(const unsigned int b) const
Definition: fe.h:3225
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
Definition: fe.cc:670
const std::vector< Point< dim > > & get_unit_support_points() const
Definition: fe.cc:1009
unsigned int global_dof_index
Definition: types.h:89
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
Definition: fe.h:3144
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
Definition: fe.cc:1260
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:940
Table< 2, int > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2516
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:3110
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const
Definition: fe.cc:972
const FiniteElement< dim, spacedim > & operator[](const unsigned int fe_index) const
Definition: fe.h:3098
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > face_system_to_base_index(const unsigned int index) const
Definition: fe.h:3213
static ::ExceptionBase & ExcEmbeddingVoid()
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:159
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
Definition: fe.h:3201
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
Definition: fe.h:3258
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2536
virtual bool hp_constraints_are_implemented() const
Definition: fe.cc:821
size_type block_start(const unsigned int i) const
Definition: fe.h:38
bool prolongation_is_implemented() const
Definition: fe.cc:696
const unsigned int first_line_index
Definition: fe_base.h:248
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:264
const bool cached_primitivity
Definition: fe.h:2634
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:205
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2618
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Definition: fe.cc:276
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2567
std::vector< std::pair< unsigned int, unsigned int > > face_system_to_component_table
Definition: fe.h:2548
virtual ~FiniteElement() override=default
bool isotropic_restriction_is_implemented() const
Definition: fe.cc:780
BlockIndices base_to_block_indices
Definition: fe.h:2579
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > face_system_to_base_table
Definition: fe.h:2573
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:253
unsigned int n_base_elements() const
Definition: fe.h:3125
unsigned int size() const
static const unsigned int space_dimension
Definition: fe.h:654
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:634
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
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:2531
TableIndices< 2 > interface_constraints_size() const
Definition: fe.cc:857
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
Definition: fe.cc:1194
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2609
UpdateFlags update_each
Definition: fe.h:715
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > operator^(const unsigned int multiplicity) const
Definition: fe.cc:150
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
Definition: fe.cc:1132
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
const std::vector< Point< dim - 1 > > & get_generalized_face_support_points() const
Definition: fe.cc:1094