Reference documentation for deal.II version Git f81eda9982 2020-03-28 21:30:57 -0400
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
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 {
3114  Assert(is_primitive(index),
3116  index)));
3117  return system_to_component_table[index];
3118 }
3119 
3120 
3121 
3122 template <int dim, int spacedim>
3123 inline unsigned int
3125 {
3126  return base_to_block_indices.size();
3127 }
3128 
3129 
3130 
3131 template <int dim, int spacedim>
3132 inline unsigned int
3134  const unsigned int index) const
3135 {
3136  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3137 }
3138 
3139 
3140 
3141 template <int dim, int spacedim>
3142 inline unsigned int
3144  const unsigned int component,
3145  const unsigned int index) const
3146 {
3147  AssertIndexRange(component, this->n_components());
3148  const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3149  std::find(system_to_component_table.begin(),
3151  std::pair<unsigned int, unsigned int>(component, index));
3152 
3153  Assert(it != system_to_component_table.end(),
3154  ExcMessage("You are asking for the number of the shape function "
3155  "within a system element that corresponds to vector "
3156  "component " +
3157  Utilities::int_to_string(component) +
3158  " and within this to "
3159  "index " +
3160  Utilities::int_to_string(index) +
3161  ". But no such "
3162  "shape function exists."));
3163  return std::distance(system_to_component_table.begin(), it);
3164 }
3165 
3166 
3167 
3168 template <int dim, int spacedim>
3169 inline std::pair<unsigned int, unsigned int>
3171  const unsigned int index) const
3172 {
3174 
3175  // in debug mode, check whether the
3176  // function is primitive, since
3177  // otherwise the result may have no
3178  // meaning
3179  //
3180  // since the primitivity tables are
3181  // all geared towards cell dof
3182  // indices, rather than face dof
3183  // indices, we have to work a
3184  // little bit...
3185  //
3186  // in 1d, the face index is equal
3187  // to the cell index
3188  Assert(is_primitive(this->face_to_cell_index(index, 0)),
3190  index)));
3191 
3192  return face_system_to_component_table[index];
3193 }
3194 
3195 
3196 
3197 template <int dim, int spacedim>
3198 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3200  const unsigned int index) const
3201 {
3202  AssertIndexRange(index, system_to_base_table.size());
3203  return system_to_base_table[index];
3204 }
3205 
3206 
3207 
3208 template <int dim, int spacedim>
3209 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3211  const unsigned int index) const
3212 {
3214  return face_system_to_base_table[index];
3215 }
3216 
3217 
3218 
3219 template <int dim, int spacedim>
3222  const unsigned int index) const
3223 {
3224  return base_to_block_indices.block_start(index);
3225 }
3226 
3227 
3228 
3229 template <int dim, int spacedim>
3230 inline std::pair<unsigned int, unsigned int>
3232  const unsigned int index) const
3233 {
3235 
3236  return component_to_base_table[index].first;
3237 }
3238 
3239 
3240 
3241 template <int dim, int spacedim>
3242 inline std::pair<unsigned int, unsigned int>
3244  const unsigned int index) const
3245 {
3246  return base_to_block_indices.global_to_local(index);
3247 }
3248 
3249 
3250 
3251 template <int dim, int spacedim>
3252 inline std::pair<unsigned int, types::global_dof_index>
3254  const unsigned int index) const
3255 {
3256  AssertIndexRange(index, this->dofs_per_cell);
3257  // The block is computed simply as
3258  // first block of this base plus
3259  // the index within the base blocks
3260  return std::pair<unsigned int, types::global_dof_index>(
3261  first_block_of_base(system_to_base_table[index].first.first) +
3262  system_to_base_table[index].first.second,
3263  system_to_base_table[index].second);
3264 }
3265 
3266 
3267 
3268 template <int dim, int spacedim>
3269 inline bool
3271  const unsigned int index) const
3272 {
3273  AssertIndexRange(index, this->dofs_per_cell);
3274  return restriction_is_additive_flags[index];
3275 }
3276 
3277 
3278 
3279 template <int dim, int spacedim>
3280 inline const ComponentMask &
3282 {
3283  AssertIndexRange(i, this->dofs_per_cell);
3284  return nonzero_components[i];
3285 }
3286 
3287 
3288 
3289 template <int dim, int spacedim>
3290 inline unsigned int
3292 {
3293  AssertIndexRange(i, this->dofs_per_cell);
3294  return n_nonzero_components_table[i];
3295 }
3296 
3297 
3298 
3299 template <int dim, int spacedim>
3300 inline bool
3302 {
3303  return cached_primitivity;
3304 }
3305 
3306 
3307 
3308 template <int dim, int spacedim>
3309 inline bool
3311 {
3312  AssertIndexRange(i, this->dofs_per_cell);
3313 
3314  // return primitivity of a shape
3315  // function by checking whether it
3316  // has more than one non-zero
3317  // component or not. we could cache
3318  // this value in an array of bools,
3319  // but accessing a bit-vector (as
3320  // std::vector<bool> is) is
3321  // probably more expensive than
3322  // just comparing against 1
3323  //
3324  // for good measure, short circuit the test
3325  // if the entire FE is primitive
3326  return (is_primitive() || (n_nonzero_components_table[i] == 1));
3327 }
3328 
3329 
3330 
3331 template <int dim, int spacedim>
3332 inline GeometryPrimitive
3334  const unsigned int cell_dof_index) const
3335 {
3336  AssertIndexRange(cell_dof_index, this->dofs_per_cell);
3337 
3338  // just go through the usual cases, taking into account how DoFs
3339  // are enumerated on the reference cell
3340  if (cell_dof_index < this->first_line_index)
3342  else if (cell_dof_index < this->first_quad_index)
3343  return GeometryPrimitive::line;
3344  else if (cell_dof_index < this->first_hex_index)
3345  return GeometryPrimitive::quad;
3346  else
3347  return GeometryPrimitive::hex;
3348 }
3349 
3350 
3351 
3352 DEAL_II_NAMESPACE_CLOSE
3353 
3354 #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:451
bool restriction_is_additive(const unsigned int index) const
Definition: fe.h:3270
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:215
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:878
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
Definition: fe.cc:967
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index) const
Definition: fe.cc:1099
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:228
FullMatrix< double > interface_constraints
Definition: fe.h:2468
bool operator!=(const FiniteElement< dim, spacedim > &) const
Definition: fe.cc:982
unsigned int n_nonzero_components(const unsigned int i) const
Definition: fe.h:3291
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:1122
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:169
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:862
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
virtual Point< dim > unit_support_point(const unsigned int index) const
Definition: fe.cc:1038
const std::vector< Point< dim - 1 > > & get_unit_face_support_points() const
Definition: fe.cc:1050
bool restriction_is_implemented() const
Definition: fe.cc:706
const std::vector< Point< dim > > & get_generalized_support_points() const
Definition: fe.cc:1016
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:181
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:239
#define AssertIndexRange(index, range)
Definition: exceptions.h:1649
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:933
bool isotropic_prolongation_is_implemented() const
Definition: fe.cc:734
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:1029
std::vector< Point< dim - 1 > > unit_face_support_points
Definition: fe.h:2487
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:365
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:1225
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:328
bool is_primitive() const
Definition: fe.h:3301
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
bool has_face_support_points() const
Definition: fe.cc:1066
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2480
bool has_support_points() const
Definition: fe.cc:1007
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:1258
unsigned int component_to_block_index(const unsigned int component) const
Definition: fe.cc:353
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:812
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:192
const unsigned int first_quad_index
Definition: fe_base.h:253
#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:3231
#define Assert(cond, exc)
Definition: exceptions.h:1419
unsigned int element_multiplicity(const unsigned int index) const
Definition: fe.h:3133
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:305
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const final
Definition: fe.cc:944
Abstract base class for mapping classes.
Definition: mapping.h:302
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
Definition: fe.h:3333
Definition: fe.h:44
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:3281
#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:894
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:790
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:3170
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:1276
virtual std::string get_name() const =0
static ::ExceptionBase & ExcProjectionVoid()
bool has_generalized_face_support_points() const
Definition: fe.cc:1090
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:911
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:537
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:431
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
Definition: fe.h:3243
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:1187
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:3221
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
Definition: fe.cc:653
const std::vector< Point< dim > > & get_unit_support_points() const
Definition: fe.cc:991
unsigned int global_dof_index
Definition: types.h:77
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
Definition: fe.h:3143
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:1240
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:922
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:954
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:3210
static ::ExceptionBase & ExcEmbeddingVoid()
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:158
unsigned int n_components() const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
Definition: fe.h:3199
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
Definition: fe.h:3253
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:803
size_type block_start(const unsigned int i) const
Definition: fe.h:38
bool prolongation_is_implemented() const
Definition: fe.cc:678
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:263
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:204
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:275
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:762
BlockIndices base_to_block_indices
Definition: fe.h:2579
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:252
unsigned int n_base_elements() const
Definition: fe.h:3124
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:618
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:839
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:1174
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:149
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
Definition: fe.cc:1112
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:1075