Reference documentation for deal.II version Git d51799cb54 2020-09-28 09:22:08 +0200
\(\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 - 2020 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_base.h>
26 #include <deal.II/fe/mapping.h>
27 
29 
30 #include <memory>
31 
32 
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 
644 template <int dim, int spacedim = dim>
645 class FiniteElement : public Subscriptor, public FiniteElementData<dim>
646 {
647 public:
651  static const unsigned int space_dimension = spacedim;
652 
678  {
679  public:
685 
689  virtual ~InternalDataBase() = default;
690 
694  InternalDataBase(const InternalDataBase &) = delete;
695 
711 
715  virtual std::size_t
716  memory_consumption() const;
717  };
718 
719 public:
762  FiniteElement(const FiniteElementData<dim> & fe_data,
763  const std::vector<bool> & restriction_is_additive_flags,
764  const std::vector<ComponentMask> &nonzero_components);
765 
769  FiniteElement(FiniteElement<dim, spacedim> &&) = default; // NOLINT
770 
774  FiniteElement(const FiniteElement<dim, spacedim> &) = default;
775 
780  virtual ~FiniteElement() override = default;
781 
790  std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>
791  operator^(const unsigned int multiplicity) const;
792 
804  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
805  clone() const = 0;
806 
817  virtual std::string
818  get_name() const = 0;
819 
842  operator[](const unsigned int fe_index) const;
843 
869  virtual double
870  shape_value(const unsigned int i, const Point<dim> &p) const;
871 
878  virtual double
879  shape_value_component(const unsigned int i,
880  const Point<dim> & p,
881  const unsigned int component) const;
882 
904  virtual Tensor<1, dim>
905  shape_grad(const unsigned int i, const Point<dim> &p) const;
906 
913  virtual Tensor<1, dim>
914  shape_grad_component(const unsigned int i,
915  const Point<dim> & p,
916  const unsigned int component) const;
917 
939  virtual Tensor<2, dim>
940  shape_grad_grad(const unsigned int i, const Point<dim> &p) const;
941 
948  virtual Tensor<2, dim>
949  shape_grad_grad_component(const unsigned int i,
950  const Point<dim> & p,
951  const unsigned int component) const;
952 
974  virtual Tensor<3, dim>
975  shape_3rd_derivative(const unsigned int i, const Point<dim> &p) const;
976 
983  virtual Tensor<3, dim>
984  shape_3rd_derivative_component(const unsigned int i,
985  const Point<dim> & p,
986  const unsigned int component) const;
987 
1009  virtual Tensor<4, dim>
1010  shape_4th_derivative(const unsigned int i, const Point<dim> &p) const;
1011 
1018  virtual Tensor<4, dim>
1019  shape_4th_derivative_component(const unsigned int i,
1020  const Point<dim> & p,
1021  const unsigned int component) const;
1032  virtual bool
1033  has_support_on_face(const unsigned int shape_index,
1034  const unsigned int face_index) const;
1035 
1037 
1057  virtual const FullMatrix<double> &
1058  get_restriction_matrix(const unsigned int child,
1059  const RefinementCase<dim> &refinement_case =
1061 
1091  virtual const FullMatrix<double> &
1092  get_prolongation_matrix(const unsigned int child,
1093  const RefinementCase<dim> &refinement_case =
1095 
1117  bool
1119 
1135  bool
1137 
1159  bool
1161 
1177  bool
1179 
1180 
1189  bool
1190  restriction_is_additive(const unsigned int index) const;
1191 
1203  const FullMatrix<double> &
1204  constraints(const ::internal::SubfaceCase<dim> &subface_case =
1206 
1222  bool
1224  const ::internal::SubfaceCase<dim> &subface_case =
1226 
1227 
1249  virtual bool
1251 
1252 
1264  virtual void
1266  FullMatrix<double> & matrix) const;
1268 
1286  virtual void
1289  const unsigned int face_no = 0) const;
1290 
1291 
1303  virtual void
1305  const unsigned int subface,
1307  const unsigned int face_no = 0) const;
1309 
1310 
1331  virtual std::vector<std::pair<unsigned int, unsigned int>>
1333 
1338  virtual std::vector<std::pair<unsigned int, unsigned int>>
1340 
1345  virtual std::vector<std::pair<unsigned int, unsigned int>>
1347  const unsigned int face_no = 0) const;
1348 
1366  const unsigned int codim = 0) const;
1367 
1369 
1400  virtual bool
1401  operator==(const FiniteElement<dim, spacedim> &fe) const;
1402 
1407  bool
1409 
1445  std::pair<unsigned int, unsigned int>
1446  system_to_component_index(const unsigned int index) const;
1447 
1457  unsigned int
1458  component_to_system_index(const unsigned int component,
1459  const unsigned int index) const;
1460 
1470  std::pair<unsigned int, unsigned int>
1471  face_system_to_component_index(const unsigned int index,
1472  const unsigned int face_no = 0) const;
1473 
1482  unsigned int
1483  adjust_quad_dof_index_for_face_orientation(const unsigned int index,
1484  const unsigned int face_no,
1485  const bool face_orientation,
1486  const bool face_flip,
1487  const bool face_rotation) const;
1488 
1543  virtual unsigned int
1544  face_to_cell_index(const unsigned int face_dof_index,
1545  const unsigned int face,
1546  const bool face_orientation = true,
1547  const bool face_flip = false,
1548  const bool face_rotation = false) const;
1549 
1557  unsigned int
1558  adjust_line_dof_index_for_line_orientation(const unsigned int index,
1559  const bool line_orientation) const;
1560 
1577  const ComponentMask &
1578  get_nonzero_components(const unsigned int i) const;
1579 
1590  unsigned int
1591  n_nonzero_components(const unsigned int i) const;
1592 
1601  bool
1602  is_primitive() const;
1603 
1614  bool
1615  is_primitive(const unsigned int i) const;
1616 
1628  unsigned int
1629  n_base_elements() const;
1630 
1635  virtual const FiniteElement<dim, spacedim> &
1636  base_element(const unsigned int index) const;
1637 
1644  unsigned int
1645  element_multiplicity(const unsigned int index) const;
1646 
1720  get_sub_fe(const ComponentMask &mask) const;
1721 
1730  virtual const FiniteElement<dim, spacedim> &
1731  get_sub_fe(const unsigned int first_component,
1732  const unsigned int n_selected_components) const;
1733 
1756  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1757  system_to_base_index(const unsigned int index) const;
1758 
1767  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1768  face_system_to_base_index(const unsigned int index,
1769  const unsigned int face_no = 0) const;
1770 
1776  first_block_of_base(const unsigned int b) const;
1777 
1790  std::pair<unsigned int, unsigned int>
1791  component_to_base_index(const unsigned int component) const;
1792 
1793 
1798  std::pair<unsigned int, unsigned int>
1799  block_to_base_index(const unsigned int block) const;
1800 
1804  std::pair<unsigned int, types::global_dof_index>
1805  system_to_block_index(const unsigned int component) const;
1806 
1810  unsigned int
1811  component_to_block_index(const unsigned int component) const;
1812 
1814 
1833  component_mask(const FEValuesExtractors::Scalar &scalar) const;
1834 
1848  component_mask(const FEValuesExtractors::Vector &vector) const;
1849 
1865  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1866 
1882  component_mask(const BlockMask &block_mask) const;
1883 
1904  BlockMask
1905  block_mask(const FEValuesExtractors::Scalar &scalar) const;
1906 
1923  BlockMask
1924  block_mask(const FEValuesExtractors::Vector &vector) const;
1925 
1943  BlockMask
1944  block_mask(const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1945 
1968  BlockMask
1969  block_mask(const ComponentMask &component_mask) const;
1970 
1986  virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
1987  get_constant_modes() const;
1988 
1990 
2026  const std::vector<Point<dim>> &
2027  get_unit_support_points() const;
2028 
2046  bool
2047  has_support_points() const;
2048 
2061  virtual Point<dim>
2062  unit_support_point(const unsigned int index) const;
2063 
2090  const std::vector<Point<dim - 1>> &
2091  get_unit_face_support_points(const unsigned int face_no = 0) const;
2092 
2101  bool
2102  has_face_support_points(const unsigned int face_no = 0) const;
2103 
2108  virtual Point<dim - 1>
2109  unit_face_support_point(const unsigned int index,
2110  const unsigned int face_no = 0) const;
2111 
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 
2269 
2278  virtual std::size_t
2279  memory_consumption() const;
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 Quadrature<dim - 1> & quadrature,
2771  FiniteElementRelatedData<dim, spacedim> &output_data) const;
2772 
2814  virtual std::unique_ptr<InternalDataBase>
2816  const UpdateFlags update_flags,
2817  const Mapping<dim, spacedim> &mapping,
2818  const Quadrature<dim - 1> & quadrature,
2820  spacedim>
2821  &output_data) const;
2822 
2903  virtual void
2905  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2906  const CellSimilarity::Similarity cell_similarity,
2907  const Quadrature<dim> & quadrature,
2908  const Mapping<dim, spacedim> & mapping,
2909  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2910  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2911  spacedim>
2912  & mapping_data,
2913  const InternalDataBase &fe_internal,
2915  spacedim>
2916  &output_data) const = 0;
2917 
2960  virtual void
2962  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2963  const unsigned int face_no,
2964  const Quadrature<dim - 1> & quadrature,
2965  const Mapping<dim, spacedim> & mapping,
2966  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2967  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2968  spacedim>
2969  & mapping_data,
2970  const InternalDataBase &fe_internal,
2972  spacedim>
2973  &output_data) const = 0;
2974 
3020  virtual void
3022  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3023  const unsigned int face_no,
3024  const unsigned int sub_no,
3025  const Quadrature<dim - 1> & quadrature,
3026  const Mapping<dim, spacedim> & mapping,
3027  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3028  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
3029  spacedim>
3030  & mapping_data,
3031  const InternalDataBase &fe_internal,
3033  spacedim>
3034  &output_data) const = 0;
3035 
3036  friend class InternalDataBase;
3037  friend class FEValuesBase<dim, spacedim>;
3038  friend class FEValues<dim, spacedim>;
3039  friend class FEFaceValues<dim, spacedim>;
3040  friend class FESubfaceValues<dim, spacedim>;
3041  friend class FESystem<dim, spacedim>;
3042 
3043  // explicitly check for sensible template arguments, but not on windows
3044  // because MSVC creates bogus warnings during normal compilation
3045 #ifndef DEAL_II_MSVC
3046  static_assert(dim <= spacedim,
3047  "The dimension <dim> of a FiniteElement must be less than or "
3048  "equal to the space dimension <spacedim> in which it lives.");
3049 #endif
3050 };
3051 
3052 
3053 //----------------------------------------------------------------------//
3054 
3055 
3056 template <int dim, int spacedim>
3058  operator[](const unsigned int fe_index) const
3059 {
3060  (void)fe_index;
3061  Assert(fe_index == 0,
3062  ExcMessage("A fe_index of zero is the only index allowed here"));
3063  return *this;
3064 }
3065 
3066 
3067 
3068 template <int dim, int spacedim>
3069 inline std::pair<unsigned int, unsigned int>
3071  const unsigned int index) const
3072 {
3074  Assert(is_primitive(index),
3076  index)));
3077  return system_to_component_table[index];
3078 }
3079 
3080 
3081 
3082 template <int dim, int spacedim>
3083 inline unsigned int
3085 {
3086  return base_to_block_indices.size();
3087 }
3088 
3089 
3090 
3091 template <int dim, int spacedim>
3092 inline unsigned int
3094  const unsigned int index) const
3095 {
3096  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3097 }
3098 
3099 
3100 
3101 template <int dim, int spacedim>
3102 inline unsigned int
3104  const unsigned int component,
3105  const unsigned int index) const
3106 {
3107  AssertIndexRange(component, this->n_components());
3108  const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3109  std::find(system_to_component_table.begin(),
3111  std::pair<unsigned int, unsigned int>(component, index));
3112 
3113  Assert(it != system_to_component_table.end(),
3114  ExcMessage("You are asking for the number of the shape function "
3115  "within a system element that corresponds to vector "
3116  "component " +
3117  Utilities::int_to_string(component) +
3118  " and within this to "
3119  "index " +
3120  Utilities::int_to_string(index) +
3121  ". But no such "
3122  "shape function exists."));
3123  return std::distance(system_to_component_table.begin(), it);
3124 }
3125 
3126 
3127 
3128 template <int dim, int spacedim>
3129 inline std::pair<unsigned int, unsigned int>
3131  const unsigned int index,
3132  const unsigned int face_no) const
3133 {
3135  index,
3136  face_system_to_component_table[this->n_unique_faces() == 1 ? 0 : face_no]
3137  .size());
3138 
3139  // in debug mode, check whether the
3140  // function is primitive, since
3141  // otherwise the result may have no
3142  // meaning
3143  //
3144  // since the primitivity tables are
3145  // all geared towards cell dof
3146  // indices, rather than face dof
3147  // indices, we have to work a
3148  // little bit...
3149  //
3150  // in 1d, the face index is equal
3151  // to the cell index
3152  Assert(is_primitive(this->face_to_cell_index(index, face_no)),
3154  index)));
3155 
3156  return face_system_to_component_table[this->n_unique_faces() == 1 ?
3157  0 :
3158  face_no][index];
3159 }
3160 
3161 
3162 
3163 template <int dim, int spacedim>
3164 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3166  const unsigned int index) const
3167 {
3168  AssertIndexRange(index, system_to_base_table.size());
3169  return system_to_base_table[index];
3170 }
3171 
3172 
3173 
3174 template <int dim, int spacedim>
3175 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3177  const unsigned int index,
3178  const unsigned int face_no) const
3179 {
3181  index,
3182  face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3183  .size());
3184  return face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3185  [index];
3186 }
3187 
3188 
3189 
3190 template <int dim, int spacedim>
3193  const unsigned int index) const
3194 {
3195  return base_to_block_indices.block_start(index);
3196 }
3197 
3198 
3199 
3200 template <int dim, int spacedim>
3201 inline std::pair<unsigned int, unsigned int>
3203  const unsigned int index) const
3204 {
3206 
3207  return component_to_base_table[index].first;
3208 }
3209 
3210 
3211 
3212 template <int dim, int spacedim>
3213 inline std::pair<unsigned int, unsigned int>
3215  const unsigned int index) const
3216 {
3217  return base_to_block_indices.global_to_local(index);
3218 }
3219 
3220 
3221 
3222 template <int dim, int spacedim>
3223 inline std::pair<unsigned int, types::global_dof_index>
3225  const unsigned int index) const
3226 {
3227  AssertIndexRange(index, this->n_dofs_per_cell());
3228  // The block is computed simply as
3229  // first block of this base plus
3230  // the index within the base blocks
3231  return std::pair<unsigned int, types::global_dof_index>(
3233  system_to_base_table[index].first.second,
3234  system_to_base_table[index].second);
3235 }
3236 
3237 
3238 
3239 template <int dim, int spacedim>
3240 inline bool
3242  const unsigned int index) const
3243 {
3244  AssertIndexRange(index, this->n_dofs_per_cell());
3245  return restriction_is_additive_flags[index];
3246 }
3247 
3248 
3249 
3250 template <int dim, int spacedim>
3251 inline const ComponentMask &
3253 {
3254  AssertIndexRange(i, this->n_dofs_per_cell());
3255  return nonzero_components[i];
3256 }
3257 
3258 
3259 
3260 template <int dim, int spacedim>
3261 inline unsigned int
3263 {
3264  AssertIndexRange(i, this->n_dofs_per_cell());
3265  return n_nonzero_components_table[i];
3266 }
3267 
3268 
3269 
3270 template <int dim, int spacedim>
3271 inline bool
3273 {
3274  return cached_primitivity;
3275 }
3276 
3277 
3278 
3279 template <int dim, int spacedim>
3280 inline bool
3282 {
3283  AssertIndexRange(i, this->n_dofs_per_cell());
3284 
3285  // return primitivity of a shape
3286  // function by checking whether it
3287  // has more than one non-zero
3288  // component or not. we could cache
3289  // this value in an array of bools,
3290  // but accessing a bit-vector (as
3291  // std::vector<bool> is) is
3292  // probably more expensive than
3293  // just comparing against 1
3294  //
3295  // for good measure, short circuit the test
3296  // if the entire FE is primitive
3297  return (is_primitive() || (n_nonzero_components_table[i] == 1));
3298 }
3299 
3300 
3301 
3302 template <int dim, int spacedim>
3303 inline GeometryPrimitive
3305  const unsigned int cell_dof_index) const
3306 {
3307  AssertIndexRange(cell_dof_index, this->n_dofs_per_cell());
3308 
3309  // just go through the usual cases, taking into account how DoFs
3310  // are enumerated on the reference cell
3311  if (cell_dof_index < this->get_first_line_index())
3313  else if (cell_dof_index < this->get_first_quad_index(0))
3314  return GeometryPrimitive::line;
3315  else if (cell_dof_index < this->get_first_hex_index())
3316  return GeometryPrimitive::quad;
3317  else
3318  return GeometryPrimitive::hex;
3319 }
3320 
3321 
3322 
3324 
3325 #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
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > face_system_to_component_table
Definition: fe.h:2507
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const =0
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:482
bool restriction_is_additive(const unsigned int index) const
Definition: fe.h:3241
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:243
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2401
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:538
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
Definition: fe.cc:660
virtual std::size_t memory_consumption() const
Definition: fe.cc:49
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
Definition: fe.cc:1029
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2452
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:256
unsigned int get_first_line_index() const
FullMatrix< double > interface_constraints
Definition: fe.h:2427
Contents is actually a matrix.
bool operator!=(const FiniteElement< dim, spacedim > &) const
Definition: fe.cc:1044
unsigned int n_nonzero_components(const unsigned int i) const
Definition: fe.h:3262
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:1168
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:197
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:931
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
virtual Point< dim > unit_support_point(const unsigned int index) const
Definition: fe.cc:1100
bool restriction_is_implemented() const
Definition: fe.cc:757
const std::vector< Point< dim > > & get_generalized_support_points() const
Definition: fe.cc:1078
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:209
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:267
bool has_face_support_points(const unsigned int face_no=0) const
Definition: fe.cc:1131
#define AssertIndexRange(index, range)
Definition: exceptions.h:1636
virtual ~InternalDataBase()=default
bool isotropic_prolongation_is_implemented() const
Definition: fe.cc:785
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2562
bool has_generalized_support_points() const
Definition: fe.cc:1091
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:396
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:1271
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:358
bool is_primitive() const
Definition: fe.h:3272
static ::ExceptionBase & ExcInterpolationNotImplemented()
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2587
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2439
bool has_support_points() const
Definition: fe.cc:1069
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:1305
unsigned int component_to_block_index(const unsigned int component) const
Definition: fe.cc:384
static ::ExceptionBase & ExcFENotPrimitive()
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index, const unsigned int face_no=0) const
Definition: fe.cc:1142
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:869
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:220
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2415
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const
Definition: fe.cc:947
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
Definition: fe.h:3202
unsigned int get_first_hex_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1411
unsigned int element_multiplicity(const unsigned int index) const
Definition: fe.h:3093
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:335
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
Definition: fe.h:3176
Abstract base class for mapping classes.
Definition: mapping.h:301
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
Definition: fe.h:3304
Definition: fe.h:44
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:3252
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
const std::vector< Point< dim - 1 > > & get_unit_face_support_points(const unsigned int face_no=0) const
Definition: fe.cc:1112
#define DeclException0(Exception0)
Definition: exceptions.h:470
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition: fe.h:2446
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:841
static ::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
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:1324
virtual std::string get_name() const =0
static ::ExceptionBase & ExcProjectionVoid()
std::vector< std::vector< Point< dim - 1 > > > generalized_face_support_points
Definition: fe.h:2458
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:982
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:568
unsigned int get_first_quad_index(const unsigned int quad_no=0) const
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:474
unsigned int n_unique_faces() const
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
Definition: fe.h:3214
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:1233
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Point< 2 > first
Definition: grid_out.cc:4335
types::global_dof_index first_block_of_base(const unsigned int b) const
Definition: fe.h:3192
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
Definition: fe.cc:704
const std::vector< Point< dim > > & get_unit_support_points() const
Definition: fe.cc:1053
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
Definition: fe.h:3103
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
Definition: fe.cc:964
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:1286
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:993
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:3070
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const
Definition: fe.cc:1016
const FiniteElement< dim, spacedim > & operator[](const unsigned int fe_index) const
Definition: fe.h:3058
static ::ExceptionBase & ExcEmbeddingVoid()
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:186
unsigned int n_components() const
unsigned int n_dofs_per_cell() const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
Definition: fe.h:3165
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
Definition: fe.h:3224
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2495
virtual bool hp_constraints_are_implemented() const
Definition: fe.cc:860
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > face_system_to_base_table
Definition: fe.h:2533
size_type block_start(const unsigned int i) const
Definition: fe.h:38
bool prolongation_is_implemented() const
Definition: fe.cc:729
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:291
const bool cached_primitivity
Definition: fe.h:2594
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:232
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2578
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Definition: fe.cc:303
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2526
virtual ~FiniteElement() override=default
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
Definition: fe.cc:1004
bool isotropic_restriction_is_implemented() const
Definition: fe.cc:813
BlockIndices base_to_block_indices
Definition: fe.h:2539
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:280
unsigned int n_base_elements() const
Definition: fe.h:3084
unsigned int size() const
static const unsigned int space_dimension
Definition: fe.h:651
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:2490
TableIndices< 2 > interface_constraints_size() const
Definition: fe.cc:903
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:1220
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2569
UpdateFlags update_each
Definition: fe.h:710
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > operator^(const unsigned int multiplicity) const
Definition: fe.cc:177
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index, const unsigned int face_no=0) const
Definition: fe.h:3130
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
Definition: fe.cc:1158
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
std::vector< Table< 2, int > > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2475