Reference documentation for deal.II version Git 6b6ce30 2018-10-20 08:13:43 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2018 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 
646 template <int dim, int spacedim = dim>
647 class FiniteElement : public Subscriptor, public FiniteElementData<dim>
648 {
649 public:
653  static const unsigned int space_dimension = spacedim;
654 
682  {
683  private:
687  InternalDataBase(const InternalDataBase &) = delete;
688 
689  public:
695 
699  virtual ~InternalDataBase() = default;
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  const FiniteElement<dim, spacedim> &fe_other) const;
1363 
1365 
1396  virtual bool
1397  operator==(const FiniteElement<dim, spacedim> &fe) const;
1398 
1403  bool
1405 
1441  std::pair<unsigned int, unsigned int>
1442  system_to_component_index(const unsigned int index) const;
1443 
1453  unsigned int
1454  component_to_system_index(const unsigned int component,
1455  const unsigned int index) const;
1456 
1466  std::pair<unsigned int, unsigned int>
1467  face_system_to_component_index(const unsigned int index) const;
1468 
1477  unsigned int
1478  adjust_quad_dof_index_for_face_orientation(const unsigned int index,
1479  const bool face_orientation,
1480  const bool face_flip,
1481  const bool face_rotation) const;
1482 
1537  virtual unsigned int
1538  face_to_cell_index(const unsigned int face_dof_index,
1539  const unsigned int face,
1540  const bool face_orientation = true,
1541  const bool face_flip = false,
1542  const bool face_rotation = false) const;
1543 
1551  unsigned int
1552  adjust_line_dof_index_for_line_orientation(const unsigned int index,
1553  const bool line_orientation) const;
1554 
1571  const ComponentMask &
1572  get_nonzero_components(const unsigned int i) const;
1573 
1584  unsigned int
1585  n_nonzero_components(const unsigned int i) const;
1586 
1595  bool
1596  is_primitive() const;
1597 
1608  bool
1609  is_primitive(const unsigned int i) const;
1610 
1622  unsigned int
1623  n_base_elements() const;
1624 
1629  virtual const FiniteElement<dim, spacedim> &
1630  base_element(const unsigned int index) const;
1631 
1638  unsigned int
1639  element_multiplicity(const unsigned int index) const;
1640 
1714  get_sub_fe(const ComponentMask &mask) const;
1715 
1724  virtual const FiniteElement<dim, spacedim> &
1725  get_sub_fe(const unsigned int first_component,
1726  const unsigned int n_selected_components) const;
1727 
1750  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1751  system_to_base_index(const unsigned int index) const;
1752 
1761  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1762  face_system_to_base_index(const unsigned int index) const;
1763 
1769  first_block_of_base(const unsigned int b) const;
1770 
1783  std::pair<unsigned int, unsigned int>
1784  component_to_base_index(const unsigned int component) const;
1785 
1786 
1791  std::pair<unsigned int, unsigned int>
1792  block_to_base_index(const unsigned int block) const;
1793 
1797  std::pair<unsigned int, types::global_dof_index>
1798  system_to_block_index(const unsigned int component) const;
1799 
1803  unsigned int
1804  component_to_block_index(const unsigned int component) const;
1805 
1807 
1826  component_mask(const FEValuesExtractors::Scalar &scalar) const;
1827 
1841  component_mask(const FEValuesExtractors::Vector &vector) const;
1842 
1858  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1859 
1875  component_mask(const BlockMask &block_mask) const;
1876 
1897  BlockMask
1898  block_mask(const FEValuesExtractors::Scalar &scalar) const;
1899 
1916  BlockMask
1917  block_mask(const FEValuesExtractors::Vector &vector) const;
1918 
1936  BlockMask
1937  block_mask(const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1938 
1961  BlockMask
1962  block_mask(const ComponentMask &component_mask) const;
1963 
1979  virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
1980  get_constant_modes() const;
1981 
1983 
2019  const std::vector<Point<dim>> &
2020  get_unit_support_points() const;
2021 
2039  bool
2040  has_support_points() const;
2041 
2054  virtual Point<dim>
2055  unit_support_point(const unsigned int index) const;
2056 
2083  const std::vector<Point<dim - 1>> &
2085 
2094  bool
2095  has_face_support_points() const;
2096 
2101  virtual Point<dim - 1>
2102  unit_face_support_point(const unsigned int index) const;
2103 
2116  const std::vector<Point<dim>> &
2118 
2128  bool
2130 
2139  DEAL_II_DEPRECATED
2140  const std::vector<Point<dim - 1>> &
2142 
2155  DEAL_II_DEPRECATED
2156  bool
2158 
2203  get_associated_geometry_primitive(const unsigned int cell_dof_index) const;
2204 
2205 
2283  virtual void
2285  const std::vector<Vector<double>> &support_point_values,
2286  std::vector<double> & nodal_values) const;
2287 
2289 
2298  virtual std::size_t
2299  memory_consumption() const;
2300 
2307  int,
2308  << "The shape function with index " << arg1
2309  << " is not primitive, i.e. it is vector-valued and "
2310  << "has more than one non-zero vector component. This "
2311  << "function cannot be called for these shape functions. "
2312  << "Maybe you want to use the same function with the "
2313  << "_component suffix?");
2327  "You are trying to access the values or derivatives of shape functions "
2328  "on the reference cell of an element that does not define its shape "
2329  "functions through mapping from the reference cell. Consequently, "
2330  "you cannot ask for shape function values or derivatives on the "
2331  "reference cell.");
2332 
2340  "You are trying to access the support points of a finite "
2341  "element that either has no support points at all, or for "
2342  "which the corresponding tables have not been implemented.");
2343 
2351  "You are trying to access the matrices that describe how "
2352  "to embed a finite element function on one cell into the "
2353  "finite element space on one of its children (i.e., the "
2354  "'embedding' or 'prolongation' matrices). However, the "
2355  "current finite element can either not define this sort of "
2356  "operation, or it has not yet been implemented.");
2357 
2366  "You are trying to access the matrices that describe how "
2367  "to restrict a finite element function from the children "
2368  "of one cell to the finite element space defined on their "
2369  "parent (i.e., the 'restriction' or 'projection' matrices). "
2370  "However, the current finite element can either not define "
2371  "this sort of operation, or it has not yet been "
2372  "implemented.");
2373 
2379  int,
2380  int,
2381  << "The interface matrix has a size of " << arg1 << "x" << arg2
2382  << ", which is not reasonable for the current element "
2383  "in the present dimension.");
2389 
2390 protected:
2404  void
2406  const bool isotropic_restriction_only = false,
2407  const bool isotropic_prolongation_only = false);
2408 
2421  std::vector<std::vector<FullMatrix<double>>> restriction;
2422 
2435  std::vector<std::vector<FullMatrix<double>>> prolongation;
2436 
2448 
2459  std::vector<Point<dim>> unit_support_points;
2460 
2466  std::vector<Point<dim - 1>> unit_face_support_points;
2467 
2472  std::vector<Point<dim>> generalized_support_points;
2473 
2479 
2496 
2511 
2515  std::vector<std::pair<unsigned int, unsigned int>> system_to_component_table;
2516 
2526  std::vector<std::pair<unsigned int, unsigned int>>
2528 
2545  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2547 
2551  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2553 
2559 
2580  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2582 
2588  const std::vector<bool> restriction_is_additive_flags;
2589 
2597  const std::vector<ComponentMask> nonzero_components;
2598 
2606  const std::vector<unsigned int> n_nonzero_components_table;
2607 
2614 
2628 
2634  static std::vector<unsigned int>
2636  const std::vector<ComponentMask> &nonzero_components);
2637 
2658  virtual UpdateFlags
2659  requires_update_flags(const UpdateFlags update_flags) const = 0;
2660 
2737  virtual std::unique_ptr<InternalDataBase>
2738  get_data(const UpdateFlags update_flags,
2739  const Mapping<dim, spacedim> &mapping,
2740  const Quadrature<dim> & quadrature,
2742  FiniteElementRelatedData<dim, spacedim> &output_data) const = 0;
2743 
2785  virtual std::unique_ptr<InternalDataBase>
2786  get_face_data(const UpdateFlags update_flags,
2787  const Mapping<dim, spacedim> &mapping,
2788  const Quadrature<dim - 1> & quadrature,
2790  FiniteElementRelatedData<dim, spacedim> &output_data) const;
2791 
2833  virtual std::unique_ptr<InternalDataBase>
2835  const UpdateFlags update_flags,
2836  const Mapping<dim, spacedim> &mapping,
2837  const Quadrature<dim - 1> & quadrature,
2839  spacedim>
2840  &output_data) const;
2841 
2922  virtual void
2924  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2925  const CellSimilarity::Similarity cell_similarity,
2926  const Quadrature<dim> & quadrature,
2927  const Mapping<dim, spacedim> & mapping,
2928  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2929  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2930  spacedim>
2931  & mapping_data,
2932  const InternalDataBase &fe_internal,
2934  spacedim>
2935  &output_data) const = 0;
2936 
2979  virtual void
2981  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2982  const unsigned int face_no,
2983  const Quadrature<dim - 1> & quadrature,
2984  const Mapping<dim, spacedim> & mapping,
2985  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2986  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2987  spacedim>
2988  & mapping_data,
2989  const InternalDataBase &fe_internal,
2991  spacedim>
2992  &output_data) const = 0;
2993 
3039  virtual void
3041  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3042  const unsigned int face_no,
3043  const unsigned int sub_no,
3044  const Quadrature<dim - 1> & quadrature,
3045  const Mapping<dim, spacedim> & mapping,
3046  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3047  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
3048  spacedim>
3049  & mapping_data,
3050  const InternalDataBase &fe_internal,
3052  spacedim>
3053  &output_data) const = 0;
3054 
3055  friend class InternalDataBase;
3056  friend class FEValuesBase<dim, spacedim>;
3057  friend class FEValues<dim, spacedim>;
3058  friend class FEFaceValues<dim, spacedim>;
3059  friend class FESubfaceValues<dim, spacedim>;
3060  friend class FESystem<dim, spacedim>;
3061 
3062  // explicitly check for sensible template arguments, but not on windows
3063  // because MSVC creates bogus warnings during normal compilation
3064 #ifndef DEAL_II_MSVC
3065  static_assert(dim <= spacedim,
3066  "The dimension <dim> of a FiniteElement must be less than or "
3067  "equal to the space dimension <spacedim> in which it lives.");
3068 #endif
3069 };
3070 
3071 
3072 //----------------------------------------------------------------------//
3073 
3074 
3075 template <int dim, int spacedim>
3077  operator[](const unsigned int fe_index) const
3078 {
3079  (void)fe_index;
3080  Assert(fe_index == 0,
3081  ExcMessage("A fe_index of zero is the only index allowed here"));
3082  return *this;
3083 }
3084 
3085 
3086 
3087 template <int dim, int spacedim>
3088 inline std::pair<unsigned int, unsigned int>
3090  const unsigned int index) const
3091 {
3092  Assert(index < system_to_component_table.size(),
3093  ExcIndexRange(index, 0, system_to_component_table.size()));
3094  Assert(is_primitive(index),
3096  index)));
3097  return system_to_component_table[index];
3098 }
3099 
3100 
3101 
3102 template <int dim, int spacedim>
3103 inline unsigned int
3105 {
3106  return base_to_block_indices.size();
3107 }
3108 
3109 
3110 
3111 template <int dim, int spacedim>
3112 inline unsigned int
3114  const unsigned int index) const
3115 {
3116  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3117 }
3118 
3119 
3120 
3121 template <int dim, int spacedim>
3122 inline unsigned int
3124  const unsigned int component,
3125  const unsigned int index) const
3126 {
3127  AssertIndexRange(component, this->n_components());
3128  const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3129  std::find(system_to_component_table.begin(),
3131  std::pair<unsigned int, unsigned int>(component, index));
3132 
3133  Assert(it != system_to_component_table.end(),
3134  ExcMessage("You are asking for the number of the shape function "
3135  "within a system element that corresponds to vector "
3136  "component " +
3137  Utilities::int_to_string(component) +
3138  " and within this to "
3139  "index " +
3140  Utilities::int_to_string(index) +
3141  ". But no such "
3142  "shape function exists."));
3143  return std::distance(system_to_component_table.begin(), it);
3144 }
3145 
3146 
3147 
3148 template <int dim, int spacedim>
3149 inline std::pair<unsigned int, unsigned int>
3151  const unsigned int index) const
3152 {
3153  Assert(index < face_system_to_component_table.size(),
3154  ExcIndexRange(index, 0, face_system_to_component_table.size()));
3155 
3156  // in debug mode, check whether the
3157  // function is primitive, since
3158  // otherwise the result may have no
3159  // meaning
3160  //
3161  // since the primitivity tables are
3162  // all geared towards cell dof
3163  // indices, rather than face dof
3164  // indices, we have to work a
3165  // little bit...
3166  //
3167  // in 1d, the face index is equal
3168  // to the cell index
3169  Assert(is_primitive(this->face_to_cell_index(index, 0)),
3171  index)));
3172 
3173  return face_system_to_component_table[index];
3174 }
3175 
3176 
3177 
3178 template <int dim, int spacedim>
3179 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3181  const unsigned int index) const
3182 {
3183  Assert(index < system_to_base_table.size(),
3184  ExcIndexRange(index, 0, system_to_base_table.size()));
3185  return system_to_base_table[index];
3186 }
3187 
3188 
3189 
3190 template <int dim, int spacedim>
3191 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3193  const unsigned int index) const
3194 {
3195  Assert(index < face_system_to_base_table.size(),
3196  ExcIndexRange(index, 0, face_system_to_base_table.size()));
3197  return face_system_to_base_table[index];
3198 }
3199 
3200 
3201 
3202 template <int dim, int spacedim>
3205  const unsigned int index) const
3206 {
3207  return base_to_block_indices.block_start(index);
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  Assert(index < component_to_base_table.size(),
3218  ExcIndexRange(index, 0, component_to_base_table.size()));
3219 
3220  return component_to_base_table[index].first;
3221 }
3222 
3223 
3224 
3225 template <int dim, int spacedim>
3226 inline std::pair<unsigned int, unsigned int>
3228  const unsigned int index) const
3229 {
3230  return base_to_block_indices.global_to_local(index);
3231 }
3232 
3233 
3234 
3235 template <int dim, int spacedim>
3236 inline std::pair<unsigned int, types::global_dof_index>
3238  const unsigned int index) const
3239 {
3240  Assert(index < this->dofs_per_cell,
3241  ExcIndexRange(index, 0, this->dofs_per_cell));
3242  // The block is computed simply as
3243  // first block of this base plus
3244  // the index within the base blocks
3245  return std::pair<unsigned int, types::global_dof_index>(
3246  first_block_of_base(system_to_base_table[index].first.first) +
3247  system_to_base_table[index].first.second,
3248  system_to_base_table[index].second);
3249 }
3250 
3251 
3252 
3253 template <int dim, int spacedim>
3254 inline bool
3256  const unsigned int index) const
3257 {
3258  Assert(index < this->dofs_per_cell,
3259  ExcIndexRange(index, 0, this->dofs_per_cell));
3260  return restriction_is_additive_flags[index];
3261 }
3262 
3263 
3264 
3265 template <int dim, int spacedim>
3266 inline const ComponentMask &
3268 {
3269  Assert(i < this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
3270  return nonzero_components[i];
3271 }
3272 
3273 
3274 
3275 template <int dim, int spacedim>
3276 inline unsigned int
3278 {
3279  Assert(i < this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
3280  return n_nonzero_components_table[i];
3281 }
3282 
3283 
3284 
3285 template <int dim, int spacedim>
3286 inline bool
3288 {
3289  return cached_primitivity;
3290 }
3291 
3292 
3293 
3294 template <int dim, int spacedim>
3295 inline bool
3297 {
3298  Assert(i < this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
3299 
3300  // return primitivity of a shape
3301  // function by checking whether it
3302  // has more than one non-zero
3303  // component or not. we could cache
3304  // this value in an array of bools,
3305  // but accessing a bit-vector (as
3306  // std::vector<bool> is) is
3307  // probably more expensive than
3308  // just comparing against 1
3309  //
3310  // for good measure, short circuit the test
3311  // if the entire FE is primitive
3312  return (is_primitive() || (n_nonzero_components_table[i] == 1));
3313 }
3314 
3315 
3316 
3317 template <int dim, int spacedim>
3318 inline GeometryPrimitive
3320  const unsigned int cell_dof_index) const
3321 {
3322  Assert(cell_dof_index < this->dofs_per_cell,
3323  ExcIndexRange(cell_dof_index, 0, this->dofs_per_cell));
3324 
3325  // just go through the usual cases, taking into account how DoFs
3326  // are enumerated on the reference cell
3327  if (cell_dof_index < this->first_line_index)
3329  else if (cell_dof_index < this->first_quad_index)
3330  return GeometryPrimitive::line;
3331  else if (cell_dof_index < this->first_hex_index)
3332  return GeometryPrimitive::quad;
3333  else
3334  return GeometryPrimitive::hex;
3335 }
3336 
3337 
3338 
3339 DEAL_II_NAMESPACE_CLOSE
3340 
3341 #endif
static::ExceptionBase & ExcFEHasNoSupportPoints()
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:159
const unsigned int first_hex_index
Definition: fe_base.h:263
bool prolongation_is_implemented() const
Definition: fe.cc:696
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const =0
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
Definition: fe.cc:1185
const std::vector< Point< dim > > & get_unit_support_points() const
Definition: fe.cc:1000
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2421
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:420
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
Definition: fe.h:3237
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:942
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2472
FullMatrix< double > interface_constraints
Definition: fe.h:2447
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:57
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
Definition: fe.h:3123
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:229
static::ExceptionBase & ExcUnitShapeValuesDoNotExist()
bool isotropic_restriction_is_implemented() const
Definition: fe.cc:780
const std::vector< Point< dim-1 > > & get_generalized_face_support_points() const
Definition: fe.cc:1085
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:882
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:808
#define AssertIndexRange(index, range)
Definition: exceptions.h:1407
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:953
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2581
std::vector< Point< dim-1 > > unit_face_support_points
Definition: fe.h:2466
virtual bool hp_constraints_are_implemented() const
Definition: fe.cc:821
bool isotropic_prolongation_is_implemented() const
Definition: fe.cc:752
const FiniteElement< dim, spacedim > & get_sub_fe(const ComponentMask &mask) const
Definition: fe.cc:1133
TableIndices< 2 > interface_constraints_size() const
Definition: fe.cc:857
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const =0
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:205
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:1236
bool operator!=(const FiniteElement< dim, spacedim > &) const
Definition: fe.cc:991
bool is_primitive() const
Definition: fe.h:3287
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
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
Definition: fe.cc:1123
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
static::ExceptionBase & ExcInterpolationNotImplemented()
std::vector< Point< dim-1 > > generalized_face_support_points
Definition: fe.h:2478
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2606
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:465
virtual std::unique_ptr< InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature,::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:931
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2459
bool restriction_is_implemented() const
Definition: fe.cc:724
unsigned int element_multiplicity(const unsigned int index) const
Definition: fe.h:3113
size_type block_start(const unsigned int i) const
static::ExceptionBase & ExcFENotPrimitive()
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
static::ExceptionBase & ExcMessage(std::string arg1)
const unsigned int first_quad_index
Definition: fe_base.h:258
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:408
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2435
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
Definition: fe.h:3227
#define Assert(cond, exc)
Definition: exceptions.h:1239
unsigned int component_to_block_index(const unsigned int component) const
Definition: fe.cc:366
UpdateFlags
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:193
unsigned int n_nonzero_components(const unsigned int i) const
Definition: fe.h:3277
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:170
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
Abstract base class for mapping classes.
Definition: dof_tools.h:57
Definition: fe.h:44
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
Definition: fe.h:3214
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:397
#define DeclException0(Exception0)
Definition: exceptions.h:385
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
Definition: fe.cc:670
static::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
Definition: fe.h:3319
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:3267
virtual Point< dim > unit_support_point(const unsigned int index) const
Definition: fe.cc:1047
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:830
virtual std::string get_name() const =0
unsigned int n_base_elements() const
Definition: fe.h:3104
static::ExceptionBase & ExcProjectionVoid()
bool has_generalized_support_points() const
Definition: fe.cc:1038
virtual Point< dim-1 > unit_face_support_point(const unsigned int index) const
Definition: fe.cc:1109
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
Definition: fe.cc:976
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index) const
Definition: fe.h:3150
virtual std::size_t memory_consumption() const
Definition: fe.cc:1215
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
Definition: fe.h:3180
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:182
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:272
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:898
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
Definition: fe.cc:914
size_type block_size(const unsigned int i) const
const FiniteElement< dim, spacedim > & operator[](const unsigned int fe_index) const
Definition: fe.h:3077
const unsigned int dofs_per_cell
Definition: fe_base.h:287
static::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
const std::vector< Point< dim > > & get_generalized_support_points() const
Definition: fe.cc:1025
unsigned int global_dof_index
Definition: types.h:89
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:1251
bool restriction_is_additive(const unsigned int index) const
Definition: fe.h:3255
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
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:1269
Table< 2, int > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2495
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
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
static::ExceptionBase & ExcEmbeddingVoid()
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > face_system_to_base_index(const unsigned int index) const
Definition: fe.h:3192
virtual std::size_t memory_consumption() const
Definition: fe.cc:49
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:3089
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:964
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2515
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:253
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:379
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
Definition: fe.h:38
const unsigned int first_line_index
Definition: fe_base.h:253
const bool cached_primitivity
Definition: fe.h:2613
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2597
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Definition: fe.cc:276
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
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2546
bool has_generalized_face_support_points() const
Definition: fe.cc:1100
unsigned int size() const
std::vector< std::pair< unsigned int, unsigned int > > face_system_to_component_table
Definition: fe.h:2527
virtual ~FiniteElement() override=default
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > operator^(const unsigned int multiplicity) const
Definition: fe.cc:150
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 std::vector< Point< dim-1 > > & get_unit_face_support_points() const
Definition: fe.cc:1060
BlockIndices base_to_block_indices
Definition: fe.h:2558
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:2552
bool has_support_points() const
Definition: fe.cc:1016
types::global_dof_index first_block_of_base(const unsigned int b) const
Definition: fe.h:3204
static const unsigned int space_dimension
Definition: fe.h:653
bool has_face_support_points() const
Definition: fe.cc:1076
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
Definition: fe.cc:1287
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:2510
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2588
UpdateFlags update_each
Definition: fe.h:715
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:1198
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