Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
fe.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1998 - 2023 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
23#include <deal.II/fe/fe_data.h>
26#include <deal.II/fe/mapping.h>
28
30#include <deal.II/lac/vector.h>
31
32#include <memory>
33
34
36
37// Forward declarations:
38#ifndef DOXYGEN
39template <int dim, int spacedim>
40class FEValuesBase;
41template <int dim, int spacedim>
42class FEValues;
43template <int dim, int spacedim>
44class FEFaceValues;
45template <int dim, int spacedim>
46class FESubfaceValues;
47namespace NonMatching
48{
49 template <int dim>
50 class FEImmersedSurfaceValues;
51}
52template <int dim, int spacedim>
53class FESystem;
54#endif
55
654template <int dim, int spacedim = dim>
655class FiniteElement : public Subscriptor, public FiniteElementData<dim>
656{
657public:
661 static constexpr unsigned int space_dimension = spacedim;
662
688 {
689 public:
695
699 virtual ~InternalDataBase() = default;
700
705
721
725 virtual std::size_t
727 };
728
729public:
773 const std::vector<bool> & restriction_is_additive_flags,
774 const std::vector<ComponentMask> &nonzero_components);
775
780
785
790 virtual ~FiniteElement() override = default;
791
800 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>
801 operator^(const unsigned int multiplicity) const;
802
814 virtual std::unique_ptr<FiniteElement<dim, spacedim>>
815 clone() const = 0;
816
827 virtual std::string
828 get_name() const = 0;
829
855 virtual double
856 shape_value(const unsigned int i, const Point<dim> &p) const;
857
864 virtual double
865 shape_value_component(const unsigned int i,
866 const Point<dim> & p,
867 const unsigned int component) const;
868
890 virtual Tensor<1, dim>
891 shape_grad(const unsigned int i, const Point<dim> &p) const;
892
899 virtual Tensor<1, dim>
900 shape_grad_component(const unsigned int i,
901 const Point<dim> & p,
902 const unsigned int component) const;
903
925 virtual Tensor<2, dim>
926 shape_grad_grad(const unsigned int i, const Point<dim> &p) const;
927
934 virtual Tensor<2, dim>
935 shape_grad_grad_component(const unsigned int i,
936 const Point<dim> & p,
937 const unsigned int component) const;
938
960 virtual Tensor<3, dim>
961 shape_3rd_derivative(const unsigned int i, const Point<dim> &p) const;
962
969 virtual Tensor<3, dim>
970 shape_3rd_derivative_component(const unsigned int i,
971 const Point<dim> & p,
972 const unsigned int component) const;
973
995 virtual Tensor<4, dim>
996 shape_4th_derivative(const unsigned int i, const Point<dim> &p) const;
997
1004 virtual Tensor<4, dim>
1006 const Point<dim> & p,
1007 const unsigned int component) const;
1018 virtual bool
1019 has_support_on_face(const unsigned int shape_index,
1020 const unsigned int face_index) const;
1021
1043 virtual const FullMatrix<double> &
1044 get_restriction_matrix(const unsigned int child,
1045 const RefinementCase<dim> &refinement_case =
1047
1077 virtual const FullMatrix<double> &
1078 get_prolongation_matrix(const unsigned int child,
1079 const RefinementCase<dim> &refinement_case =
1081
1103 bool
1105
1121 bool
1123
1145 bool
1147
1163 bool
1165
1166
1175 bool
1176 restriction_is_additive(const unsigned int index) const;
1177
1189 const FullMatrix<double> &
1190 constraints(const ::internal::SubfaceCase<dim> &subface_case =
1192
1208 bool
1210 const ::internal::SubfaceCase<dim> &subface_case =
1212
1213
1235 virtual bool
1237
1238
1250 virtual void
1252 FullMatrix<double> & matrix) const;
1272 virtual void
1274 FullMatrix<double> & matrix,
1275 const unsigned int face_no = 0) const;
1276
1277
1289 virtual void
1291 const unsigned int subface,
1292 FullMatrix<double> & matrix,
1293 const unsigned int face_no = 0) const;
1317 virtual std::vector<std::pair<unsigned int, unsigned int>>
1319
1324 virtual std::vector<std::pair<unsigned int, unsigned int>>
1326
1331 virtual std::vector<std::pair<unsigned int, unsigned int>>
1333 const unsigned int face_no = 0) const;
1334
1352 const unsigned int codim = 0) const;
1353
1386 virtual bool
1388
1393 bool
1395
1431 std::pair<unsigned int, unsigned int>
1432 system_to_component_index(const unsigned int index) const;
1433
1443 unsigned int
1444 component_to_system_index(const unsigned int component,
1445 const unsigned int index) const;
1446
1456 std::pair<unsigned int, unsigned int>
1457 face_system_to_component_index(const unsigned int index,
1458 const unsigned int face_no = 0) const;
1459
1468 unsigned int
1470 const unsigned int face_no,
1471 const bool face_orientation,
1472 const bool face_flip,
1473 const bool face_rotation) const;
1474
1529 virtual unsigned int
1530 face_to_cell_index(const unsigned int face_dof_index,
1531 const unsigned int face,
1532 const bool face_orientation = true,
1533 const bool face_flip = false,
1534 const bool face_rotation = false) const;
1535
1543 unsigned int
1545 const bool line_orientation) const;
1546
1563 const ComponentMask &
1564 get_nonzero_components(const unsigned int i) const;
1565
1576 unsigned int
1577 n_nonzero_components(const unsigned int i) const;
1578
1587 bool
1589
1600 bool
1601 is_primitive(const unsigned int i) const;
1602
1614 unsigned int
1616
1621 virtual const FiniteElement<dim, spacedim> &
1622 base_element(const unsigned int index) const;
1623
1630 unsigned int
1631 element_multiplicity(const unsigned int index) const;
1632
1706 get_sub_fe(const ComponentMask &mask) const;
1707
1716 virtual const FiniteElement<dim, spacedim> &
1717 get_sub_fe(const unsigned int first_component,
1718 const unsigned int n_selected_components) const;
1719
1742 std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1743 system_to_base_index(const unsigned int index) const;
1744
1753 std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1754 face_system_to_base_index(const unsigned int index,
1755 const unsigned int face_no = 0) const;
1756
1762 first_block_of_base(const unsigned int b) const;
1763
1776 std::pair<unsigned int, unsigned int>
1777 component_to_base_index(const unsigned int component) const;
1778
1779
1784 std::pair<unsigned int, unsigned int>
1785 block_to_base_index(const unsigned int block) const;
1786
1790 std::pair<unsigned int, types::global_dof_index>
1791 system_to_block_index(const unsigned int component) const;
1792
1796 unsigned int
1797 component_to_block_index(const unsigned int component) const;
1798
1820
1835
1851 const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1852
1869
1890 BlockMask
1892
1909 BlockMask
1911
1929 BlockMask
1931
1954 BlockMask
1956
1972 virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
1974
2012 const std::vector<Point<dim>> &
2014
2032 bool
2034
2047 virtual Point<dim>
2048 unit_support_point(const unsigned int index) const;
2049
2076 const std::vector<Point<dim - 1>> &
2077 get_unit_face_support_points(const unsigned int face_no = 0) const;
2078
2087 bool
2088 has_face_support_points(const unsigned int face_no = 0) const;
2089
2094 virtual Point<dim - 1>
2095 unit_face_support_point(const unsigned int index,
2096 const unsigned int face_no = 0) const;
2097
2125 const std::vector<Point<dim>> &
2127
2137 bool
2139
2184 get_associated_geometry_primitive(const unsigned int cell_dof_index) const;
2185
2186
2264 virtual void
2266 const std::vector<Vector<double>> &support_point_values,
2267 std::vector<double> & nodal_values) const;
2268
2279 virtual std::size_t
2281
2288 int,
2289 << "The shape function with index " << arg1
2290 << " is not primitive, i.e. it is vector-valued and "
2291 << "has more than one non-zero vector component. This "
2292 << "function cannot be called for these shape functions. "
2293 << "Maybe you want to use the same function with the "
2294 << "_component suffix?");
2308 "You are trying to access the values or derivatives of shape functions "
2309 "on the reference cell of an element that does not define its shape "
2310 "functions through mapping from the reference cell. Consequently, "
2311 "you cannot ask for shape function values or derivatives on the "
2312 "reference cell.");
2313
2321 "You are trying to access the support points of a finite "
2322 "element that either has no support points at all, or for "
2323 "which the corresponding tables have not been implemented.");
2324
2332 "You are trying to access the matrices that describe how "
2333 "to embed a finite element function on one cell into the "
2334 "finite element space on one of its children (i.e., the "
2335 "'embedding' or 'prolongation' matrices). However, the "
2336 "current finite element can either not define this sort of "
2337 "operation, or it has not yet been implemented.");
2338
2347 "You are trying to access the matrices that describe how "
2348 "to restrict a finite element function from the children "
2349 "of one cell to the finite element space defined on their "
2350 "parent (i.e., the 'restriction' or 'projection' matrices). "
2351 "However, the current finite element can either not define "
2352 "this sort of operation, or it has not yet been "
2353 "implemented.");
2354
2360 int,
2361 int,
2362 << "The interface matrix has a size of " << arg1 << 'x' << arg2
2363 << ", which is not reasonable for the current element "
2364 "in the present dimension.");
2370
2371protected:
2385 void
2387 const bool isotropic_restriction_only = false,
2388 const bool isotropic_prolongation_only = false);
2389
2402 std::vector<std::vector<FullMatrix<double>>> restriction;
2403
2416 std::vector<std::vector<FullMatrix<double>>> prolongation;
2417
2429
2440 std::vector<Point<dim>> unit_support_points;
2441
2447 std::vector<std::vector<Point<dim - 1>>> unit_face_support_points;
2448
2453 std::vector<Point<dim>> generalized_support_points;
2454
2459 std::vector<std::vector<Point<dim - 1>>> generalized_face_support_points;
2460
2477
2492
2496 std::vector<std::pair<unsigned int, unsigned int>> system_to_component_table;
2497
2507 std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
2509
2526 std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2528
2532 std::vector<
2533 std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>>
2535
2541
2562 std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2564
2570 const std::vector<bool> restriction_is_additive_flags;
2571
2579 const std::vector<ComponentMask> nonzero_components;
2580
2588 const std::vector<unsigned int> n_nonzero_components_table;
2589
2596
2610
2616 static std::vector<unsigned int>
2618 const std::vector<ComponentMask> &nonzero_components);
2619
2640 virtual UpdateFlags
2641 requires_update_flags(const UpdateFlags update_flags) const = 0;
2642
2719 virtual std::unique_ptr<InternalDataBase>
2720 get_data(const UpdateFlags update_flags,
2721 const Mapping<dim, spacedim> &mapping,
2722 const Quadrature<dim> & quadrature,
2723 ::internal::FEValuesImplementation::
2724 FiniteElementRelatedData<dim, spacedim> &output_data) const = 0;
2725
2767 virtual std::unique_ptr<InternalDataBase>
2768 get_face_data(const UpdateFlags update_flags,
2769 const Mapping<dim, spacedim> & mapping,
2770 const hp::QCollection<dim - 1> &quadrature,
2771 ::internal::FEValuesImplementation::
2772 FiniteElementRelatedData<dim, spacedim> &output_data) const;
2773
2777 virtual std::unique_ptr<InternalDataBase>
2779 const UpdateFlags update_flags,
2780 const Mapping<dim, spacedim> &mapping,
2781 const Quadrature<dim - 1> & quadrature,
2783 &output_data) const;
2784
2826 virtual std::unique_ptr<InternalDataBase>
2828 const UpdateFlags update_flags,
2829 const Mapping<dim, spacedim> &mapping,
2830 const Quadrature<dim - 1> & quadrature,
2832 spacedim>
2833 &output_data) const;
2834
2915 virtual void
2918 const CellSimilarity::Similarity cell_similarity,
2919 const Quadrature<dim> & quadrature,
2920 const Mapping<dim, spacedim> & mapping,
2921 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2923 & mapping_data,
2924 const InternalDataBase &fe_internal,
2926 spacedim>
2927 &output_data) const = 0;
2928
2971 virtual void
2974 const unsigned int face_no,
2975 const hp::QCollection<dim - 1> & quadrature,
2976 const Mapping<dim, spacedim> & mapping,
2977 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2979 & mapping_data,
2980 const InternalDataBase &fe_internal,
2982 spacedim>
2983 &output_data) const;
2984
2988 virtual void
2991 const unsigned int face_no,
2992 const Quadrature<dim - 1> & quadrature,
2993 const Mapping<dim, spacedim> & mapping,
2994 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
2996 & mapping_data,
2997 const InternalDataBase &fe_internal,
2999 &output_data) const;
3000
3046 virtual void
3049 const unsigned int face_no,
3050 const unsigned int sub_no,
3051 const Quadrature<dim - 1> & quadrature,
3052 const Mapping<dim, spacedim> & mapping,
3053 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3055 & mapping_data,
3056 const InternalDataBase &fe_internal,
3058 spacedim>
3059 &output_data) const = 0;
3060
3061 friend class InternalDataBase;
3062 friend class FEValuesBase<dim, spacedim>;
3063 friend class FEValues<dim, spacedim>;
3064 friend class FEFaceValues<dim, spacedim>;
3065 friend class FESubfaceValues<dim, spacedim>;
3066 friend class NonMatching::FEImmersedSurfaceValues<dim>;
3067 friend class FESystem<dim, spacedim>;
3068
3069 // explicitly check for sensible template arguments, but not on windows
3070 // because MSVC creates bogus warnings during normal compilation
3071#ifndef DEAL_II_MSVC
3072 static_assert(dim <= spacedim,
3073 "The dimension <dim> of a FiniteElement must be less than or "
3074 "equal to the space dimension <spacedim> in which it lives.");
3075#endif
3076};
3077
3078
3079//----------------------------------------------------------------------//
3080
3081
3082template <int dim, int spacedim>
3083inline std::pair<unsigned int, unsigned int>
3085 const unsigned int index) const
3086{
3088 Assert(is_primitive(index),
3090 index)));
3091 return system_to_component_table[index];
3092}
3093
3094
3095
3096template <int dim, int spacedim>
3097inline unsigned int
3099{
3100 return base_to_block_indices.size();
3101}
3102
3103
3104
3105template <int dim, int spacedim>
3106inline unsigned int
3108 const unsigned int index) const
3109{
3110 return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3111}
3112
3113
3114
3115template <int dim, int spacedim>
3116inline unsigned int
3118 const unsigned int component,
3119 const unsigned int index) const
3120{
3121 AssertIndexRange(component, this->n_components());
3122 const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3123 std::find(system_to_component_table.begin(),
3125 std::pair<unsigned int, unsigned int>(component, index));
3126
3128 ExcMessage("You are asking for the number of the shape function "
3129 "within a system element that corresponds to vector "
3130 "component " +
3131 Utilities::int_to_string(component) +
3132 " and within this to "
3133 "index " +
3135 ". But no such "
3136 "shape function exists."));
3137 return std::distance(system_to_component_table.begin(), it);
3138}
3139
3140
3141
3142template <int dim, int spacedim>
3143inline std::pair<unsigned int, unsigned int>
3145 const unsigned int index,
3146 const unsigned int face_no) const
3147{
3149 index,
3150 face_system_to_component_table[this->n_unique_faces() == 1 ? 0 : face_no]
3151 .size());
3152
3153 // in debug mode, check whether the
3154 // function is primitive, since
3155 // otherwise the result may have no
3156 // meaning
3157 //
3158 // since the primitivity tables are
3159 // all geared towards cell dof
3160 // indices, rather than face dof
3161 // indices, we have to work a
3162 // little bit...
3163 //
3164 // in 1d, the face index is equal
3165 // to the cell index
3166 Assert(is_primitive(this->face_to_cell_index(index, face_no)),
3168 index)));
3169
3171 0 :
3172 face_no][index];
3173}
3174
3175
3176
3177template <int dim, int spacedim>
3178inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3180 const unsigned int index) const
3181{
3184}
3185
3186
3187
3188template <int dim, int spacedim>
3189inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3191 const unsigned int index,
3192 const unsigned int face_no) const
3193{
3195 index,
3196 face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3197 .size());
3198 return face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3199 [index];
3200}
3201
3202
3203
3204template <int dim, int spacedim>
3207 const unsigned int index) const
3208{
3209 return base_to_block_indices.block_start(index);
3210}
3211
3212
3213
3214template <int dim, int spacedim>
3215inline std::pair<unsigned int, unsigned int>
3217 const unsigned int index) const
3218{
3220
3221 return component_to_base_table[index].first;
3222}
3223
3224
3225
3226template <int dim, int spacedim>
3227inline std::pair<unsigned int, unsigned int>
3229 const unsigned int index) const
3230{
3232}
3233
3234
3235
3236template <int dim, int spacedim>
3237inline std::pair<unsigned int, types::global_dof_index>
3239 const unsigned int index) const
3240{
3241 AssertIndexRange(index, this->n_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>(
3247 system_to_base_table[index].first.second,
3248 system_to_base_table[index].second);
3249}
3250
3251
3252
3253template <int dim, int spacedim>
3254inline bool
3256 const unsigned int index) const
3257{
3258 AssertIndexRange(index, this->n_dofs_per_cell());
3259 return restriction_is_additive_flags[index];
3260}
3261
3262
3263
3264template <int dim, int spacedim>
3265inline const ComponentMask &
3267{
3269 return nonzero_components[i];
3270}
3271
3272
3273
3274template <int dim, int spacedim>
3275inline unsigned int
3277{
3279 return n_nonzero_components_table[i];
3280}
3281
3282
3283
3284template <int dim, int spacedim>
3285inline bool
3287{
3288 return cached_primitivity;
3289}
3290
3291
3292
3293template <int dim, int spacedim>
3294inline bool
3295FiniteElement<dim, spacedim>::is_primitive(const unsigned int i) const
3296{
3298
3299 // return primitivity of a shape
3300 // function by checking whether it
3301 // has more than one non-zero
3302 // component or not. we could cache
3303 // this value in an array of bools,
3304 // but accessing a bit-vector (as
3305 // std::vector<bool> is) is
3306 // probably more expensive than
3307 // just comparing against 1
3308 //
3309 // for good measure, short circuit the test
3310 // if the entire FE is primitive
3311 return (is_primitive() || (n_nonzero_components_table[i] == 1));
3312}
3313
3314
3315
3316template <int dim, int spacedim>
3317inline GeometryPrimitive
3319 const unsigned int cell_dof_index) const
3320{
3321 AssertIndexRange(cell_dof_index, this->n_dofs_per_cell());
3322
3323 // just go through the usual cases, taking into account how DoFs
3324 // are enumerated on the reference cell
3325 if (cell_dof_index < this->get_first_line_index())
3327 else if (cell_dof_index < this->get_first_quad_index(0))
3329 else if (cell_dof_index < this->get_first_hex_index())
3331 else
3333}
3334
3335
3336
3338
3339#endif
size_type block_size(const unsigned int i) const
unsigned int size() const
size_type block_start(const unsigned int i) const
std::pair< unsigned int, size_type > global_to_local(const size_type i) const
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
unsigned int get_first_line_index() const
unsigned int n_dofs_per_cell() const
unsigned int get_first_quad_index(const unsigned int quad_no=0) const
unsigned int n_components() const
unsigned int n_unique_faces() const
unsigned int get_first_hex_index() const
virtual ~InternalDataBase()=default
InternalDataBase(const InternalDataBase &)=delete
virtual std::size_t memory_consumption() const
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
bool isotropic_prolongation_is_implemented() const
virtual std::string get_name() const =0
virtual Point< dim > unit_support_point(const unsigned int index) const
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const
const std::vector< unsigned int > n_nonzero_components_table
Definition fe.h:2588
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 const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
const ComponentMask & get_nonzero_components(const unsigned int i) const
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< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition fe.h:2447
bool prolongation_is_implemented() const
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
ComponentMask component_mask(const BlockMask &block_mask) const
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
ComponentMask component_mask(const FEValuesExtractors::SymmetricTensor< 2 > &sym_tensor) const
virtual std::unique_ptr< InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
std::pair< unsigned int, unsigned int > component_to_base_index(const unsigned int component) const
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index, const unsigned int face_no=0) const
bool is_primitive(const unsigned int i) const
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
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
bool has_face_support_points(const unsigned int face_no=0) const
const bool cached_primitivity
Definition fe.h:2595
virtual const FiniteElement< dim, spacedim > & get_sub_fe(const unsigned int first_component, const unsigned int n_selected_components) const
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
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
bool operator!=(const FiniteElement< dim, spacedim > &) const
ComponentMask component_mask(const FEValuesExtractors::Vector &vector) const
bool has_support_points() const
std::vector< std::vector< FullMatrix< double > > > restriction
Definition fe.h:2402
bool has_generalized_support_points() const
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
std::vector< Table< 2, int > > adjust_quad_dof_index_for_face_orientation_table
Definition fe.h:2476
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
const std::vector< Point< dim > > & get_unit_support_points() const
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< 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
static constexpr unsigned int space_dimension
Definition fe.h:661
BlockIndices base_to_block_indices
Definition fe.h:2540
bool is_primitive() const
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
const std::vector< Point< dim - 1 > > & get_unit_face_support_points(const unsigned int face_no=0) const
bool isotropic_restriction_is_implemented() const
std::vector< std::vector< Point< dim - 1 > > > generalized_face_support_points
Definition fe.h:2459
std::pair< unsigned int, unsigned int > block_to_base_index(const unsigned int block) const
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition fe.h:2491
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
bool restriction_is_additive(const unsigned int index) const
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
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 FiniteElement< dim, spacedim > & get_sub_fe(const ComponentMask &mask) const
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
FiniteElement(const FiniteElement< dim, spacedim > &)=default
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
unsigned int component_to_block_index(const unsigned int component) const
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const =0
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition fe.h:2527
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double > > &support_point_values, std::vector< double > &nodal_values) const
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
BlockMask block_mask(const FEValuesExtractors::Vector &vector) const
BlockMask block_mask(const ComponentMask &component_mask) const
std::pair< std::pair< unsigned int, unsigned int >, unsigned int > system_to_base_index(const unsigned int index) const
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition fe.h:2496
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
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
unsigned int element_multiplicity(const unsigned int index) const
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > operator^(const unsigned int multiplicity) const
const std::vector< bool > restriction_is_additive_flags
Definition fe.h:2570
virtual ~FiniteElement() override=default
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition fe.h:2563
TableIndices< 2 > interface_constraints_size() const
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > face_system_to_base_table
Definition fe.h:2534
std::vector< Point< dim > > unit_support_points
Definition fe.h:2440
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
bool restriction_is_implemented() const
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
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > face_system_to_component_table
Definition fe.h:2508
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
FiniteElement(FiniteElement< dim, spacedim > &&)=default
unsigned int n_nonzero_components(const unsigned int i) const
FullMatrix< double > interface_constraints
Definition fe.h:2428
const std::vector< Point< dim > > & get_generalized_support_points() const
BlockMask block_mask(const FEValuesExtractors::SymmetricTensor< 2 > &sym_tensor) const
unsigned int n_base_elements() const
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const
virtual std::size_t memory_consumption() const
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
GeometryPrimitive get_associated_geometry_primitive(const unsigned int cell_dof_index) const
std::vector< Point< dim > > generalized_support_points
Definition fe.h:2453
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
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
const std::vector< ComponentMask > nonzero_components
Definition fe.h:2579
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index, const unsigned int face_no=0) const
virtual bool hp_constraints_are_implemented() const
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition fe.h:2416
types::global_dof_index first_block_of_base(const unsigned int b) const
Abstract base class for mapping classes.
Definition mapping.h:317
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 2 > first
Definition grid_out.cc:4615
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:533
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
static ::ExceptionBase & ExcProjectionVoid()
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcEmbeddingVoid()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
static ::ExceptionBase & ExcInterpolationNotImplemented()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcFEHasNoSupportPoints()
UpdateFlags
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition utilities.cc:471