Reference documentation for deal.II version Git b43ef918fe 2022-01-29 10:49:23 +0100
\(\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 // Forward declarations:
36 #ifndef DOXYGEN
37 template <int dim, int spacedim>
38 class FEValuesBase;
39 template <int dim, int spacedim>
40 class FEValues;
41 template <int dim, int spacedim>
42 class FEFaceValues;
43 template <int dim, int spacedim>
44 class FESubfaceValues;
45 namespace NonMatching
46 {
47  template <int dim>
48  class FEImmersedSurfaceValues;
49 }
50 template <int dim, int spacedim>
51 class FESystem;
52 #endif
53 
652 template <int dim, int spacedim = dim>
653 class FiniteElement : public Subscriptor, public FiniteElementData<dim>
654 {
655 public:
659  static const unsigned int space_dimension = spacedim;
660 
686  {
687  public:
693 
697  virtual ~InternalDataBase() = default;
698 
702  InternalDataBase(const InternalDataBase &) = delete;
703 
719 
723  virtual std::size_t
724  memory_consumption() const;
725  };
726 
727 public:
770  FiniteElement(const FiniteElementData<dim> & fe_data,
771  const std::vector<bool> & restriction_is_additive_flags,
772  const std::vector<ComponentMask> &nonzero_components);
773 
777  FiniteElement(FiniteElement<dim, spacedim> &&) = default; // NOLINT
778 
782  FiniteElement(const FiniteElement<dim, spacedim> &) = default;
783 
788  virtual ~FiniteElement() override = default;
789 
798  std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>
799  operator^(const unsigned int multiplicity) const;
800 
812  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
813  clone() const = 0;
814 
825  virtual std::string
826  get_name() const = 0;
827 
853  operator[](const unsigned int fe_index) const;
854 
880  virtual double
881  shape_value(const unsigned int i, const Point<dim> &p) const;
882 
889  virtual double
890  shape_value_component(const unsigned int i,
891  const Point<dim> & p,
892  const unsigned int component) const;
893 
915  virtual Tensor<1, dim>
916  shape_grad(const unsigned int i, const Point<dim> &p) const;
917 
924  virtual Tensor<1, dim>
925  shape_grad_component(const unsigned int i,
926  const Point<dim> & p,
927  const unsigned int component) const;
928 
950  virtual Tensor<2, dim>
951  shape_grad_grad(const unsigned int i, const Point<dim> &p) const;
952 
959  virtual Tensor<2, dim>
960  shape_grad_grad_component(const unsigned int i,
961  const Point<dim> & p,
962  const unsigned int component) const;
963 
985  virtual Tensor<3, dim>
986  shape_3rd_derivative(const unsigned int i, const Point<dim> &p) const;
987 
994  virtual Tensor<3, dim>
995  shape_3rd_derivative_component(const unsigned int i,
996  const Point<dim> & p,
997  const unsigned int component) const;
998 
1020  virtual Tensor<4, dim>
1021  shape_4th_derivative(const unsigned int i, const Point<dim> &p) const;
1022 
1029  virtual Tensor<4, dim>
1030  shape_4th_derivative_component(const unsigned int i,
1031  const Point<dim> & p,
1032  const unsigned int component) const;
1043  virtual bool
1044  has_support_on_face(const unsigned int shape_index,
1045  const unsigned int face_index) const;
1046 
1048 
1068  virtual const FullMatrix<double> &
1069  get_restriction_matrix(const unsigned int child,
1070  const RefinementCase<dim> &refinement_case =
1072 
1102  virtual const FullMatrix<double> &
1103  get_prolongation_matrix(const unsigned int child,
1104  const RefinementCase<dim> &refinement_case =
1106 
1128  bool
1129  prolongation_is_implemented() const;
1130 
1146  bool
1147  isotropic_prolongation_is_implemented() const;
1148 
1170  bool
1171  restriction_is_implemented() const;
1172 
1188  bool
1189  isotropic_restriction_is_implemented() const;
1190 
1191 
1200  bool
1201  restriction_is_additive(const unsigned int index) const;
1202 
1214  const FullMatrix<double> &
1215  constraints(const ::internal::SubfaceCase<dim> &subface_case =
1217 
1233  bool
1234  constraints_are_implemented(
1235  const ::internal::SubfaceCase<dim> &subface_case =
1237 
1238 
1260  virtual bool
1261  hp_constraints_are_implemented() const;
1262 
1263 
1275  virtual void
1277  FullMatrix<double> & matrix) const;
1279 
1297  virtual void
1298  get_face_interpolation_matrix(const FiniteElement<dim, spacedim> &source,
1300  const unsigned int face_no = 0) const;
1301 
1302 
1314  virtual void
1315  get_subface_interpolation_matrix(const FiniteElement<dim, spacedim> &source,
1316  const unsigned int subface,
1318  const unsigned int face_no = 0) const;
1320 
1321 
1342  virtual std::vector<std::pair<unsigned int, unsigned int>>
1343  hp_vertex_dof_identities(const FiniteElement<dim, spacedim> &fe_other) const;
1344 
1349  virtual std::vector<std::pair<unsigned int, unsigned int>>
1350  hp_line_dof_identities(const FiniteElement<dim, spacedim> &fe_other) const;
1351 
1356  virtual std::vector<std::pair<unsigned int, unsigned int>>
1357  hp_quad_dof_identities(const FiniteElement<dim, spacedim> &fe_other,
1358  const unsigned int face_no = 0) const;
1359 
1376  compare_for_domination(const FiniteElement<dim, spacedim> &fe_other,
1377  const unsigned int codim = 0) const;
1378 
1380 
1411  virtual bool
1412  operator==(const FiniteElement<dim, spacedim> &fe) const;
1413 
1418  bool
1420 
1456  std::pair<unsigned int, unsigned int>
1457  system_to_component_index(const unsigned int index) const;
1458 
1468  unsigned int
1469  component_to_system_index(const unsigned int component,
1470  const unsigned int index) const;
1471 
1481  std::pair<unsigned int, unsigned int>
1482  face_system_to_component_index(const unsigned int index,
1483  const unsigned int face_no = 0) const;
1484 
1493  unsigned int
1494  adjust_quad_dof_index_for_face_orientation(const unsigned int index,
1495  const unsigned int face_no,
1496  const bool face_orientation,
1497  const bool face_flip,
1498  const bool face_rotation) const;
1499 
1554  virtual unsigned int
1555  face_to_cell_index(const unsigned int face_dof_index,
1556  const unsigned int face,
1557  const bool face_orientation = true,
1558  const bool face_flip = false,
1559  const bool face_rotation = false) const;
1560 
1568  unsigned int
1569  adjust_line_dof_index_for_line_orientation(const unsigned int index,
1570  const bool line_orientation) const;
1571 
1588  const ComponentMask &
1589  get_nonzero_components(const unsigned int i) const;
1590 
1601  unsigned int
1602  n_nonzero_components(const unsigned int i) const;
1603 
1612  bool
1613  is_primitive() const;
1614 
1625  bool
1626  is_primitive(const unsigned int i) const;
1627 
1639  unsigned int
1640  n_base_elements() const;
1641 
1646  virtual const FiniteElement<dim, spacedim> &
1647  base_element(const unsigned int index) const;
1648 
1655  unsigned int
1656  element_multiplicity(const unsigned int index) const;
1657 
1731  get_sub_fe(const ComponentMask &mask) const;
1732 
1741  virtual const FiniteElement<dim, spacedim> &
1742  get_sub_fe(const unsigned int first_component,
1743  const unsigned int n_selected_components) const;
1744 
1767  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1768  system_to_base_index(const unsigned int index) const;
1769 
1778  std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1779  face_system_to_base_index(const unsigned int index,
1780  const unsigned int face_no = 0) const;
1781 
1787  first_block_of_base(const unsigned int b) const;
1788 
1801  std::pair<unsigned int, unsigned int>
1802  component_to_base_index(const unsigned int component) const;
1803 
1804 
1809  std::pair<unsigned int, unsigned int>
1810  block_to_base_index(const unsigned int block) const;
1811 
1815  std::pair<unsigned int, types::global_dof_index>
1816  system_to_block_index(const unsigned int component) const;
1817 
1821  unsigned int
1822  component_to_block_index(const unsigned int component) const;
1823 
1825 
1844  component_mask(const FEValuesExtractors::Scalar &scalar) const;
1845 
1859  component_mask(const FEValuesExtractors::Vector &vector) const;
1860 
1875  component_mask(
1876  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1877 
1893  component_mask(const BlockMask &block_mask) const;
1894 
1915  BlockMask
1916  block_mask(const FEValuesExtractors::Scalar &scalar) const;
1917 
1934  BlockMask
1935  block_mask(const FEValuesExtractors::Vector &vector) const;
1936 
1954  BlockMask
1955  block_mask(const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1956 
1979  BlockMask
1980  block_mask(const ComponentMask &component_mask) const;
1981 
1997  virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
1998  get_constant_modes() const;
1999 
2001 
2037  const std::vector<Point<dim>> &
2038  get_unit_support_points() const;
2039 
2057  bool
2058  has_support_points() const;
2059 
2072  virtual Point<dim>
2073  unit_support_point(const unsigned int index) const;
2074 
2101  const std::vector<Point<dim - 1>> &
2102  get_unit_face_support_points(const unsigned int face_no = 0) const;
2103 
2112  bool
2113  has_face_support_points(const unsigned int face_no = 0) const;
2114 
2119  virtual Point<dim - 1>
2120  unit_face_support_point(const unsigned int index,
2121  const unsigned int face_no = 0) const;
2122 
2135  const std::vector<Point<dim>> &
2136  get_generalized_support_points() const;
2137 
2147  bool
2148  has_generalized_support_points() const;
2149 
2194  get_associated_geometry_primitive(const unsigned int cell_dof_index) const;
2195 
2196 
2274  virtual void
2276  const std::vector<Vector<double>> &support_point_values,
2277  std::vector<double> & nodal_values) const;
2278 
2280 
2289  virtual std::size_t
2290  memory_consumption() const;
2291 
2297  DeclException1(ExcShapeFunctionNotPrimitive,
2298  int,
2299  << "The shape function with index " << arg1
2300  << " is not primitive, i.e. it is vector-valued and "
2301  << "has more than one non-zero vector component. This "
2302  << "function cannot be called for these shape functions. "
2303  << "Maybe you want to use the same function with the "
2304  << "_component suffix?");
2317  ExcUnitShapeValuesDoNotExist,
2318  "You are trying to access the values or derivatives of shape functions "
2319  "on the reference cell of an element that does not define its shape "
2320  "functions through mapping from the reference cell. Consequently, "
2321  "you cannot ask for shape function values or derivatives on the "
2322  "reference cell.");
2323 
2330  DeclExceptionMsg(ExcFEHasNoSupportPoints,
2331  "You are trying to access the support points of a finite "
2332  "element that either has no support points at all, or for "
2333  "which the corresponding tables have not been implemented.");
2334 
2341  DeclExceptionMsg(ExcEmbeddingVoid,
2342  "You are trying to access the matrices that describe how "
2343  "to embed a finite element function on one cell into the "
2344  "finite element space on one of its children (i.e., the "
2345  "'embedding' or 'prolongation' matrices). However, the "
2346  "current finite element can either not define this sort of "
2347  "operation, or it has not yet been implemented.");
2348 
2356  DeclExceptionMsg(ExcProjectionVoid,
2357  "You are trying to access the matrices that describe how "
2358  "to restrict a finite element function from the children "
2359  "of one cell to the finite element space defined on their "
2360  "parent (i.e., the 'restriction' or 'projection' matrices). "
2361  "However, the current finite element can either not define "
2362  "this sort of operation, or it has not yet been "
2363  "implemented.");
2364 
2369  DeclException2(ExcWrongInterfaceMatrixSize,
2370  int,
2371  int,
2372  << "The interface matrix has a size of " << arg1 << "x" << arg2
2373  << ", which is not reasonable for the current element "
2374  "in the present dimension.");
2379  DeclException0(ExcInterpolationNotImplemented);
2380 
2381 protected:
2395  void
2396  reinit_restriction_and_prolongation_matrices(
2397  const bool isotropic_restriction_only = false,
2398  const bool isotropic_prolongation_only = false);
2399 
2412  std::vector<std::vector<FullMatrix<double>>> restriction;
2413 
2426  std::vector<std::vector<FullMatrix<double>>> prolongation;
2427 
2439 
2450  std::vector<Point<dim>> unit_support_points;
2451 
2457  std::vector<std::vector<Point<dim - 1>>> unit_face_support_points;
2458 
2463  std::vector<Point<dim>> generalized_support_points;
2464 
2469  std::vector<std::vector<Point<dim - 1>>> generalized_face_support_points;
2470 
2487 
2502 
2506  std::vector<std::pair<unsigned int, unsigned int>> system_to_component_table;
2507 
2517  std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
2519 
2536  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2538 
2542  std::vector<
2543  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>>
2545 
2551 
2572  std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2574 
2580  const std::vector<bool> restriction_is_additive_flags;
2581 
2589  const std::vector<ComponentMask> nonzero_components;
2590 
2598  const std::vector<unsigned int> n_nonzero_components_table;
2599 
2606 
2619  interface_constraints_size() const;
2620 
2626  static std::vector<unsigned int>
2627  compute_n_nonzero_components(
2628  const std::vector<ComponentMask> &nonzero_components);
2629 
2650  virtual UpdateFlags
2651  requires_update_flags(const UpdateFlags update_flags) const = 0;
2652 
2729  virtual std::unique_ptr<InternalDataBase>
2730  get_data(const UpdateFlags update_flags,
2731  const Mapping<dim, spacedim> &mapping,
2732  const Quadrature<dim> & quadrature,
2734  FiniteElementRelatedData<dim, spacedim> &output_data) const = 0;
2735 
2777  virtual std::unique_ptr<InternalDataBase>
2778  get_face_data(const UpdateFlags update_flags,
2779  const Mapping<dim, spacedim> & mapping,
2780  const hp::QCollection<dim - 1> &quadrature,
2782  FiniteElementRelatedData<dim, spacedim> &output_data) const;
2783 
2787  virtual std::unique_ptr<InternalDataBase>
2788  get_face_data(
2789  const UpdateFlags update_flags,
2790  const Mapping<dim, spacedim> &mapping,
2791  const Quadrature<dim - 1> & quadrature,
2793  &output_data) const;
2794 
2836  virtual std::unique_ptr<InternalDataBase>
2837  get_subface_data(
2838  const UpdateFlags update_flags,
2839  const Mapping<dim, spacedim> &mapping,
2840  const Quadrature<dim - 1> & quadrature,
2842  spacedim>
2843  &output_data) const;
2844 
2925  virtual void
2926  fill_fe_values(
2927  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2928  const CellSimilarity::Similarity cell_similarity,
2929  const Quadrature<dim> & quadrature,
2930  const Mapping<dim, spacedim> & mapping,
2931  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2932  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2933  spacedim>
2934  & mapping_data,
2935  const InternalDataBase &fe_internal,
2937  spacedim>
2938  &output_data) const = 0;
2939 
2982  virtual void
2983  fill_fe_face_values(
2984  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
2985  const unsigned int face_no,
2986  const hp::QCollection<dim - 1> & quadrature,
2987  const Mapping<dim, spacedim> & mapping,
2988  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2989  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
2990  spacedim>
2991  & mapping_data,
2992  const InternalDataBase &fe_internal,
2994  spacedim>
2995  &output_data) const;
2996 
3000  virtual void
3001  fill_fe_face_values(
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,
3008  & mapping_data,
3009  const InternalDataBase &fe_internal,
3011  &output_data) const;
3012 
3058  virtual void
3059  fill_fe_subface_values(
3060  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
3061  const unsigned int face_no,
3062  const unsigned int sub_no,
3063  const Quadrature<dim - 1> & quadrature,
3064  const Mapping<dim, spacedim> & mapping,
3065  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3066  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
3067  spacedim>
3068  & mapping_data,
3069  const InternalDataBase &fe_internal,
3071  spacedim>
3072  &output_data) const = 0;
3073 
3074  friend class InternalDataBase;
3075  friend class FEValuesBase<dim, spacedim>;
3076  friend class FEValues<dim, spacedim>;
3077  friend class FEFaceValues<dim, spacedim>;
3078  friend class FESubfaceValues<dim, spacedim>;
3080  friend class FESystem<dim, spacedim>;
3081 
3082  // explicitly check for sensible template arguments, but not on windows
3083  // because MSVC creates bogus warnings during normal compilation
3084 #ifndef DEAL_II_MSVC
3085  static_assert(dim <= spacedim,
3086  "The dimension <dim> of a FiniteElement must be less than or "
3087  "equal to the space dimension <spacedim> in which it lives.");
3088 #endif
3089 };
3090 
3091 
3092 //----------------------------------------------------------------------//
3093 
3094 
3095 template <int dim, int spacedim>
3096 inline const FiniteElement<dim, spacedim> &
3097 FiniteElement<dim, spacedim>::operator[](const unsigned int fe_index) const
3098 {
3099  (void)fe_index;
3100  Assert(fe_index == 0,
3101  ExcMessage("A fe_index of zero is the only index allowed here"));
3102  return *this;
3103 }
3104 
3105 
3106 
3107 template <int dim, int spacedim>
3108 inline std::pair<unsigned int, unsigned int>
3110  const unsigned int index) const
3111 {
3113  Assert(is_primitive(index),
3115  index)));
3116  return system_to_component_table[index];
3117 }
3118 
3119 
3120 
3121 template <int dim, int spacedim>
3122 inline unsigned int
3124 {
3125  return base_to_block_indices.size();
3126 }
3127 
3128 
3129 
3130 template <int dim, int spacedim>
3131 inline unsigned int
3133  const unsigned int index) const
3134 {
3135  return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3136 }
3137 
3138 
3139 
3140 template <int dim, int spacedim>
3141 inline unsigned int
3143  const unsigned int component,
3144  const unsigned int index) const
3145 {
3146  AssertIndexRange(component, this->n_components());
3147  const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3148  std::find(system_to_component_table.begin(),
3150  std::pair<unsigned int, unsigned int>(component, index));
3151 
3152  Assert(it != system_to_component_table.end(),
3153  ExcMessage("You are asking for the number of the shape function "
3154  "within a system element that corresponds to vector "
3155  "component " +
3156  Utilities::int_to_string(component) +
3157  " and within this to "
3158  "index " +
3159  Utilities::int_to_string(index) +
3160  ". But no such "
3161  "shape function exists."));
3162  return std::distance(system_to_component_table.begin(), it);
3163 }
3164 
3165 
3166 
3167 template <int dim, int spacedim>
3168 inline std::pair<unsigned int, unsigned int>
3170  const unsigned int index,
3171  const unsigned int face_no) const
3172 {
3174  index,
3175  face_system_to_component_table[this->n_unique_faces() == 1 ? 0 : face_no]
3176  .size());
3177 
3178  // in debug mode, check whether the
3179  // function is primitive, since
3180  // otherwise the result may have no
3181  // meaning
3182  //
3183  // since the primitivity tables are
3184  // all geared towards cell dof
3185  // indices, rather than face dof
3186  // indices, we have to work a
3187  // little bit...
3188  //
3189  // in 1d, the face index is equal
3190  // to the cell index
3191  Assert(is_primitive(this->face_to_cell_index(index, face_no)),
3193  index)));
3194 
3195  return face_system_to_component_table[this->n_unique_faces() == 1 ?
3196  0 :
3197  face_no][index];
3198 }
3199 
3200 
3201 
3202 template <int dim, int spacedim>
3203 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3205  const unsigned int index) const
3206 {
3207  AssertIndexRange(index, system_to_base_table.size());
3208  return system_to_base_table[index];
3209 }
3210 
3211 
3212 
3213 template <int dim, int spacedim>
3214 inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3216  const unsigned int index,
3217  const unsigned int face_no) const
3218 {
3220  index,
3221  face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3222  .size());
3223  return face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3224  [index];
3225 }
3226 
3227 
3228 
3229 template <int dim, int spacedim>
3232  const unsigned int index) const
3233 {
3234  return base_to_block_indices.block_start(index);
3235 }
3236 
3237 
3238 
3239 template <int dim, int spacedim>
3240 inline std::pair<unsigned int, unsigned int>
3242  const unsigned int index) const
3243 {
3245 
3246  return component_to_base_table[index].first;
3247 }
3248 
3249 
3250 
3251 template <int dim, int spacedim>
3252 inline std::pair<unsigned int, unsigned int>
3254  const unsigned int index) const
3255 {
3256  return base_to_block_indices.global_to_local(index);
3257 }
3258 
3259 
3260 
3261 template <int dim, int spacedim>
3262 inline std::pair<unsigned int, types::global_dof_index>
3264  const unsigned int index) const
3265 {
3266  AssertIndexRange(index, this->n_dofs_per_cell());
3267  // The block is computed simply as
3268  // first block of this base plus
3269  // the index within the base blocks
3270  return std::pair<unsigned int, types::global_dof_index>(
3272  system_to_base_table[index].first.second,
3273  system_to_base_table[index].second);
3274 }
3275 
3276 
3277 
3278 template <int dim, int spacedim>
3279 inline bool
3281  const unsigned int index) const
3282 {
3283  AssertIndexRange(index, this->n_dofs_per_cell());
3284  return restriction_is_additive_flags[index];
3285 }
3286 
3287 
3288 
3289 template <int dim, int spacedim>
3290 inline const ComponentMask &
3292 {
3293  AssertIndexRange(i, this->n_dofs_per_cell());
3294  return nonzero_components[i];
3295 }
3296 
3297 
3298 
3299 template <int dim, int spacedim>
3300 inline unsigned int
3302 {
3303  AssertIndexRange(i, this->n_dofs_per_cell());
3304  return n_nonzero_components_table[i];
3305 }
3306 
3307 
3308 
3309 template <int dim, int spacedim>
3310 inline bool
3312 {
3313  return cached_primitivity;
3314 }
3315 
3316 
3317 
3318 template <int dim, int spacedim>
3319 inline bool
3321 {
3322  AssertIndexRange(i, this->n_dofs_per_cell());
3323 
3324  // return primitivity of a shape
3325  // function by checking whether it
3326  // has more than one non-zero
3327  // component or not. we could cache
3328  // this value in an array of bools,
3329  // but accessing a bit-vector (as
3330  // std::vector<bool> is) is
3331  // probably more expensive than
3332  // just comparing against 1
3333  //
3334  // for good measure, short circuit the test
3335  // if the entire FE is primitive
3336  return (is_primitive() || (n_nonzero_components_table[i] == 1));
3337 }
3338 
3339 
3340 
3341 template <int dim, int spacedim>
3342 inline GeometryPrimitive
3344  const unsigned int cell_dof_index) const
3345 {
3346  AssertIndexRange(cell_dof_index, this->n_dofs_per_cell());
3347 
3348  // just go through the usual cases, taking into account how DoFs
3349  // are enumerated on the reference cell
3350  if (cell_dof_index < this->get_first_line_index())
3352  else if (cell_dof_index < this->get_first_quad_index(0))
3353  return GeometryPrimitive::line;
3354  else if (cell_dof_index < this->get_first_hex_index())
3355  return GeometryPrimitive::quad;
3356  else
3357  return GeometryPrimitive::hex;
3358 }
3359 
3360 
3361 
3363 
3364 #endif
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > face_system_to_component_table
Definition: fe.h:2518
bool restriction_is_additive(const unsigned int index) const
Definition: fe.h:3280
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2412
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2463
unsigned int get_first_line_index() const
FullMatrix< double > interface_constraints
Definition: fe.h:2438
Contents is actually a matrix.
unsigned int n_nonzero_components(const unsigned int i) const
Definition: fe.h:3301
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1720
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2573
size_type block_size(const unsigned int i) const
bool is_primitive() const
Definition: fe.h:3311
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:855
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2598
void convert_generalized_support_point_values_to_dof_values(const FiniteElement< dim, spacedim > &finite_element, const std::vector< Vector< number >> &support_point_values, std::vector< number > &dof_values)
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
static ::ExceptionBase & ExcFENotPrimitive()
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2450
static ::ExceptionBase & ExcMessage(std::string arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
void get_interpolation_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &interpolation_matrix)
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2426
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
Definition: fe.h:3241
unsigned int get_first_hex_index() const
#define Assert(cond, exc)
Definition: exceptions.h:1461
unsigned int element_multiplicity(const unsigned int index) const
Definition: fe.h:3132
UpdateFlags
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:3215
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
Definition: fe.h:3343
Abstract base class for mapping classes.
Definition: mapping.h:310
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:3291
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
#define DeclException0(Exception0)
Definition: exceptions.h:464
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition: fe.h:2457
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:407
std::vector< std::vector< Point< dim - 1 > > > generalized_face_support_points
Definition: fe.h:2469
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:473
unsigned int n_unique_faces() const
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
Definition: fe.h:3253
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Point< 2 > first
Definition: grid_out.cc:4614
types::global_dof_index first_block_of_base(const unsigned int b) const
Definition: fe.h:3231
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
Definition: fe.h:3142
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:3109
const FiniteElement< dim, spacedim > & operator[](const unsigned int fe_index) const
Definition: fe.h:3097
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:3204
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
Definition: fe.h:3263
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:406
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2506
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > face_system_to_base_table
Definition: fe.h:2544
size_type block_start(const unsigned int i) const
Expression operator^(const Expression &lhs, const Expression &rhs)
const bool cached_primitivity
Definition: fe.h:2605
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2589
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2537
BlockIndices base_to_block_indices
Definition: fe.h:2550
unsigned int n_base_elements() const
Definition: fe.h:3123
unsigned int size() const
#define DEAL_II_DEPRECATED
Definition: config.h:163
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:2501
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2580
UpdateFlags update_each
Definition: fe.h:718
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:3169
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
std::vector< Table< 2, int > > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2486