deal.II version GIT relicensing-2846-g6fb608615f 2025-03-15 04:10:00+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1998 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_fe_h
16#define dealii_fe_h
17
18#include <deal.II/base/config.h>
19
22#include <deal.II/fe/fe_data.h>
25#include <deal.II/fe/mapping.h>
27
29#include <deal.II/lac/vector.h>
30
31#include <memory>
32
33
35
36// Forward declarations:
37#ifndef DOXYGEN
38template <int dim, int spacedim>
39class FEValuesBase;
40template <int dim, int spacedim>
41class FEValues;
42template <int dim, int spacedim>
43class FEFaceValues;
44template <int dim, int spacedim>
45class FESubfaceValues;
46namespace NonMatching
47{
48 template <int dim>
49 class FEImmersedSurfaceValues;
50}
51template <int dim, int spacedim>
52class FESystem;
53#endif
54
653template <int dim, int spacedim = dim>
655 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
970 virtual Tensor<3, dim>
971 shape_3rd_derivative_component(const unsigned int i,
972 const Point<dim> &p,
973 const unsigned int component) const;
974
996 virtual Tensor<4, dim>
997 shape_4th_derivative(const unsigned int i, const Point<dim> &p) const;
998
1006 virtual Tensor<4, dim>
1008 const Point<dim> &p,
1009 const unsigned int component) const;
1020 virtual bool
1021 has_support_on_face(const unsigned int shape_index,
1022 const unsigned int face_index) const;
1023
1060 virtual const FullMatrix<double> &
1061 get_restriction_matrix(const unsigned int child,
1062 const RefinementCase<dim> &refinement_case =
1064
1094 virtual const FullMatrix<double> &
1095 get_prolongation_matrix(const unsigned int child,
1096 const RefinementCase<dim> &refinement_case =
1098
1120 bool
1122
1138 bool
1140
1162 bool
1164
1180 bool
1182
1183
1192 bool
1193 restriction_is_additive(const unsigned int index) const;
1194
1206 const FullMatrix<double> &
1207 constraints(const ::internal::SubfaceCase<dim> &subface_case =
1209
1225 bool
1227 const ::internal::SubfaceCase<dim> &subface_case =
1229
1230
1252 virtual bool
1254
1255
1267 virtual void
1269 FullMatrix<double> &matrix) const;
1289 virtual void
1291 FullMatrix<double> &matrix,
1292 const unsigned int face_no = 0) const;
1293
1294
1306 virtual void
1308 const unsigned int subface,
1309 FullMatrix<double> &matrix,
1310 const unsigned int face_no = 0) const;
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 const unsigned int face_no = 0) const;
1351
1369 const unsigned int codim = 0) const;
1370
1403 virtual bool
1405
1410 bool
1412
1448 std::pair<unsigned int, unsigned int>
1449 system_to_component_index(const unsigned int index) const;
1450
1460 unsigned int
1461 component_to_system_index(const unsigned int component,
1462 const unsigned int index) const;
1463
1473 std::pair<unsigned int, unsigned int>
1474 face_system_to_component_index(const unsigned int index,
1475 const unsigned int face_no = 0) const;
1476
1483 unsigned int
1485 const unsigned int index,
1486 const unsigned int face_no,
1487 const types::geometric_orientation combined_orientation) const;
1488
1536 virtual unsigned int
1537 face_to_cell_index(const unsigned int face_dof_index,
1538 const unsigned int face,
1539 const types::geometric_orientation combined_orientation =
1541
1551 unsigned int
1553 const unsigned int index,
1554 const types::geometric_orientation combined_orientation) const;
1555
1572 const ComponentMask &
1573 get_nonzero_components(const unsigned int i) const;
1574
1585 unsigned int
1586 n_nonzero_components(const unsigned int i) const;
1587
1596 bool
1598
1609 bool
1610 is_primitive(const unsigned int i) const;
1611
1628 virtual const Table<2, bool> &
1630
1642 unsigned int
1644
1649 virtual const FiniteElement<dim, spacedim> &
1650 base_element(const unsigned int index) const;
1651
1658 unsigned int
1659 element_multiplicity(const unsigned int index) const;
1660
1734 get_sub_fe(const ComponentMask &mask) const;
1735
1744 virtual const FiniteElement<dim, spacedim> &
1745 get_sub_fe(const unsigned int first_component,
1746 const unsigned int n_selected_components) const;
1747
1770 std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1771 system_to_base_index(const unsigned int index) const;
1772
1781 std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
1782 face_system_to_base_index(const unsigned int index,
1783 const unsigned int face_no = 0) const;
1784
1790 first_block_of_base(const unsigned int b) const;
1791
1804 std::pair<unsigned int, unsigned int>
1805 component_to_base_index(const unsigned int component) const;
1806
1807
1812 std::pair<unsigned int, unsigned int>
1813 block_to_base_index(const unsigned int block) const;
1814
1818 std::pair<unsigned int, types::global_dof_index>
1819 system_to_block_index(const unsigned int component) const;
1820
1824 unsigned int
1825 component_to_block_index(const unsigned int component) const;
1826
1842 bool
1843 shape_function_belongs_to(const unsigned int shape_function,
1844 const FEValuesExtractors::Scalar &component) const;
1845
1861 bool
1862 shape_function_belongs_to(const unsigned int shape_function,
1864
1865
1881 bool
1883 const unsigned int shape_function,
1885
1886
1897 bool
1899 const unsigned int shape_function,
1922
1937
1953 const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
1954
1971
1992 BlockMask
1994
2011 BlockMask
2013
2031 BlockMask
2033
2056 BlockMask
2058
2074 virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
2076
2114 const std::vector<Point<dim>> &
2116
2143 bool
2145
2158 virtual Point<dim>
2159 unit_support_point(const unsigned int index) const;
2160
2187 const std::vector<Point<dim - 1>> &
2188 get_unit_face_support_points(const unsigned int face_no = 0) const;
2189
2201 bool
2202 has_face_support_points(const unsigned int face_no = 0) const;
2203
2208 virtual Point<dim - 1>
2209 unit_face_support_point(const unsigned int index,
2210 const unsigned int face_no = 0) const;
2211
2239 const std::vector<Point<dim>> &
2241
2256 bool
2258
2305 get_associated_geometry_primitive(const unsigned int cell_dof_index) const;
2306
2307
2385 virtual void
2387 const std::vector<Vector<double>> &support_point_values,
2388 std::vector<double> &nodal_values) const;
2389
2400 virtual std::size_t
2402
2409 int,
2410 << "The shape function with index " << arg1
2411 << " is not primitive, i.e. it is vector-valued and "
2412 << "has more than one non-zero vector component. This "
2413 << "function cannot be called for these shape functions. "
2414 << "Maybe you want to use the same function with the "
2415 << "_component suffix?");
2429 "You are trying to access the values or derivatives of shape functions "
2430 "on the reference cell of an element that does not define its shape "
2431 "functions through mapping from the reference cell. Consequently, "
2432 "you cannot ask for shape function values or derivatives on the "
2433 "reference cell.");
2434
2442 "You are trying to access the support points of a finite "
2443 "element that either has no support points at all, or for "
2444 "which the corresponding tables have not been implemented.");
2445
2453 "You are trying to access the matrices that describe how "
2454 "to embed a finite element function on one cell into the "
2455 "finite element space on one of its children (i.e., the "
2456 "'embedding' or 'prolongation' matrices). However, the "
2457 "current finite element can either not define this sort of "
2458 "operation, or it has not yet been implemented.");
2459
2468 "You are trying to access the matrices that describe how "
2469 "to restrict a finite element function from the children "
2470 "of one cell to the finite element space defined on their "
2471 "parent (i.e., the 'restriction' or 'projection' matrices). "
2472 "However, the current finite element can either not define "
2473 "this sort of operation, or it has not yet been "
2474 "implemented.");
2475
2481 int,
2482 int,
2483 << "The interface matrix has a size of " << arg1 << 'x' << arg2
2484 << ", which is not reasonable for the current element "
2485 "in the present dimension.");
2491
2492protected:
2506 void
2508 const bool isotropic_restriction_only = false,
2509 const bool isotropic_prolongation_only = false);
2510
2523 std::vector<std::vector<FullMatrix<double>>> restriction;
2524
2537 std::vector<std::vector<FullMatrix<double>>> prolongation;
2538
2550
2561 std::vector<Point<dim>> unit_support_points;
2562
2568 std::vector<std::vector<Point<dim - 1>>> unit_face_support_points;
2569
2574 std::vector<Point<dim>> generalized_support_points;
2575
2580 std::vector<std::vector<Point<dim - 1>>> generalized_face_support_points;
2581
2598
2611
2615 std::vector<std::pair<unsigned int, unsigned int>> system_to_component_table;
2616
2626 std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
2628
2645 std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2647
2651 std::vector<
2652 std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>>
2654
2660
2681 std::vector<std::pair<std::pair<unsigned int, unsigned int>, unsigned int>>
2683
2689 const std::vector<bool> restriction_is_additive_flags;
2690
2698 const std::vector<ComponentMask> nonzero_components;
2699
2707 const std::vector<unsigned int> n_nonzero_components_table;
2708
2715
2723
2738
2744 static std::vector<unsigned int>
2746 const std::vector<ComponentMask> &nonzero_components);
2747
2768 virtual UpdateFlags
2769 requires_update_flags(const UpdateFlags update_flags) const = 0;
2770
2848 virtual std::unique_ptr<InternalDataBase>
2849 get_data(const UpdateFlags update_flags,
2850 const Mapping<dim, spacedim> &mapping,
2851 const Quadrature<dim> &quadrature,
2852 ::internal::FEValuesImplementation::
2853 FiniteElementRelatedData<dim, spacedim> &output_data) const = 0;
2854
2897 virtual std::unique_ptr<InternalDataBase>
2898 get_face_data(const UpdateFlags update_flags,
2899 const Mapping<dim, spacedim> &mapping,
2900 const hp::QCollection<dim - 1> &quadrature,
2901 ::internal::FEValuesImplementation::
2902 FiniteElementRelatedData<dim, spacedim> &output_data) const;
2903
2907 virtual std::unique_ptr<InternalDataBase>
2909 const UpdateFlags update_flags,
2910 const Mapping<dim, spacedim> &mapping,
2911 const Quadrature<dim - 1> &quadrature,
2913 &output_data) const;
2914
2958 virtual std::unique_ptr<InternalDataBase>
2960 const UpdateFlags update_flags,
2961 const Mapping<dim, spacedim> &mapping,
2962 const Quadrature<dim - 1> &quadrature,
2964 spacedim>
2965 &output_data) const;
2966
3047 virtual void
3050 const CellSimilarity::Similarity cell_similarity,
3051 const Quadrature<dim> &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
3103 virtual void
3106 const unsigned int face_no,
3107 const hp::QCollection<dim - 1> &quadrature,
3108 const Mapping<dim, spacedim> &mapping,
3109 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3111 &mapping_data,
3112 const InternalDataBase &fe_internal,
3114 spacedim>
3115 &output_data) const;
3116
3120 virtual void
3123 const unsigned int face_no,
3124 const Quadrature<dim - 1> &quadrature,
3125 const Mapping<dim, spacedim> &mapping,
3126 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3128 &mapping_data,
3129 const InternalDataBase &fe_internal,
3131 &output_data) const;
3132
3178 virtual void
3181 const unsigned int face_no,
3182 const unsigned int sub_no,
3183 const Quadrature<dim - 1> &quadrature,
3184 const Mapping<dim, spacedim> &mapping,
3185 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
3187 &mapping_data,
3188 const InternalDataBase &fe_internal,
3190 spacedim>
3191 &output_data) const = 0;
3192
3193 friend class InternalDataBase;
3194 friend class FEValuesBase<dim, spacedim>;
3195 friend class FEValues<dim, spacedim>;
3196 friend class FEFaceValues<dim, spacedim>;
3197 friend class FESubfaceValues<dim, spacedim>;
3198 friend class NonMatching::FEImmersedSurfaceValues<dim>;
3199 friend class FESystem<dim, spacedim>;
3200
3201 // explicitly check for sensible template arguments, but not on windows
3202 // because MSVC creates bogus warnings during normal compilation
3203#ifndef DEAL_II_MSVC
3204 static_assert(dim <= spacedim,
3205 "The dimension <dim> of a FiniteElement must be less than or "
3206 "equal to the space dimension <spacedim> in which it lives.");
3207#endif
3208};
3209
3210
3211//----------------------------------------------------------------------//
3212#ifndef DOXYGEN
3213
3214template <int dim, int spacedim>
3215inline std::pair<unsigned int, unsigned int>
3217 const unsigned int index) const
3218{
3220 Assert(is_primitive(index),
3222 index)));
3223 return system_to_component_table[index];
3224}
3225
3226
3227
3228template <int dim, int spacedim>
3229inline unsigned int
3231{
3232 return base_to_block_indices.size();
3233}
3234
3235
3236
3237template <int dim, int spacedim>
3238inline unsigned int
3240 const unsigned int index) const
3241{
3242 return static_cast<unsigned int>(base_to_block_indices.block_size(index));
3243}
3244
3245
3246
3247template <int dim, int spacedim>
3248inline unsigned int
3250 const unsigned int component,
3251 const unsigned int index) const
3252{
3253 AssertIndexRange(component, this->n_components());
3254 const std::vector<std::pair<unsigned int, unsigned int>>::const_iterator it =
3255 std::find(system_to_component_table.begin(),
3257 std::pair<unsigned int, unsigned int>(component, index));
3258
3260 ExcMessage("You are asking for the number of the shape function "
3261 "within a system element that corresponds to vector "
3262 "component " +
3263 Utilities::int_to_string(component) +
3264 " and within this to "
3265 "index " +
3267 ". But no such "
3268 "shape function exists."));
3269 return std::distance(system_to_component_table.begin(), it);
3270}
3271
3272
3273
3274template <int dim, int spacedim>
3275inline std::pair<unsigned int, unsigned int>
3277 const unsigned int index,
3278 const unsigned int face_no) const
3279{
3281 index,
3282 face_system_to_component_table[this->n_unique_faces() == 1 ? 0 : face_no]
3283 .size());
3284
3285 // in debug mode, check whether the
3286 // function is primitive, since
3287 // otherwise the result may have no
3288 // meaning
3289 //
3290 // since the primitivity tables are
3291 // all geared towards cell dof
3292 // indices, rather than face dof
3293 // indices, we have to work a
3294 // little bit...
3295 //
3296 // in 1d, the face index is equal
3297 // to the cell index
3298 Assert(is_primitive(this->face_to_cell_index(index, face_no)),
3300 index)));
3301
3303 0 :
3304 face_no][index];
3305}
3306
3307
3308
3309template <int dim, int spacedim>
3310inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3312 const unsigned int index) const
3313{
3316}
3317
3318
3319
3320template <int dim, int spacedim>
3321inline std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
3323 const unsigned int index,
3324 const unsigned int face_no) const
3325{
3327 index,
3328 face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3329 .size());
3330 return face_system_to_base_table[this->n_unique_faces() == 1 ? 0 : face_no]
3331 [index];
3332}
3333
3334
3335
3336template <int dim, int spacedim>
3339 const unsigned int index) const
3340{
3341 return base_to_block_indices.block_start(index);
3342}
3343
3344
3345
3346template <int dim, int spacedim>
3347inline std::pair<unsigned int, unsigned int>
3349 const unsigned int index) const
3350{
3352
3353 return component_to_base_table[index].first;
3354}
3355
3356
3357
3358template <int dim, int spacedim>
3359inline std::pair<unsigned int, unsigned int>
3361 const unsigned int index) const
3362{
3364}
3365
3366
3367
3368template <int dim, int spacedim>
3369inline std::pair<unsigned int, types::global_dof_index>
3371 const unsigned int index) const
3372{
3373 AssertIndexRange(index, this->n_dofs_per_cell());
3374 // The block is computed simply as
3375 // first block of this base plus
3376 // the index within the base blocks
3377 return std::pair<unsigned int, types::global_dof_index>(
3379 system_to_base_table[index].first.second,
3380 system_to_base_table[index].second);
3381}
3382
3383
3384
3385template <int dim, int spacedim>
3386inline bool
3388 const unsigned int shape_function,
3389 const FEValuesExtractors::Scalar &component) const
3390{
3391 AssertIndexRange(shape_function, this->n_dofs_per_cell());
3392 AssertIndexRange(component.component, this->n_components());
3393
3394 if (is_primitive(shape_function))
3395 return (system_to_component_table[shape_function].first ==
3396 component.component);
3397 else
3398 return nonzero_components[shape_function][component.component];
3399}
3400
3401
3402
3403template <int dim, int spacedim>
3404inline bool
3406 const unsigned int shape_function,
3408{
3409 AssertIndexRange(shape_function, this->n_dofs_per_cell());
3410 AssertIndexRange(components.first_vector_component, this->n_components());
3411 AssertIndexRange(components.first_vector_component + dim,
3412 this->n_components() + 1);
3413
3414 if (is_primitive(shape_function))
3415 return ((system_to_component_table[shape_function].first >=
3416 components.first_vector_component) &&
3417 (system_to_component_table[shape_function].first <
3418 components.first_vector_component + dim));
3419 else
3420 // Return whether there is overlap between the nonzero components
3421 // of the current shape function and the selected components of
3422 // the extractor:
3423 {
3424 for (unsigned int i = components.first_vector_component;
3425 i < components.first_vector_component + dim;
3426 ++i)
3427 if (nonzero_components[shape_function][i])
3428 return true;
3429 return false;
3430 }
3431}
3432
3433
3434
3435template <int dim, int spacedim>
3436inline bool
3438 const unsigned int shape_function,
3440{
3441 AssertIndexRange(shape_function, this->n_dofs_per_cell());
3442 AssertIndexRange(components.first_tensor_component, this->n_components());
3443 AssertIndexRange(components.first_tensor_component + dim,
3444 this->n_components() + 1);
3445
3446 if (is_primitive(shape_function))
3447 return ((system_to_component_table[shape_function].first >=
3448 components.first_tensor_component) &&
3449 (system_to_component_table[shape_function].first <
3450 components.first_tensor_component +
3452 else
3453 // Return whether there is overlap between the nonzero components
3454 // of the current shape function and the selected components of
3455 // the extractor:
3456 {
3457 for (unsigned int i = components.first_tensor_component;
3458 i < components.first_tensor_component +
3460 ++i)
3461 if (nonzero_components[shape_function][i])
3462 return true;
3463 return false;
3464 }
3465}
3466
3467
3468
3469template <int dim, int spacedim>
3470inline bool
3472 const unsigned int shape_function,
3474{
3475 AssertIndexRange(shape_function, this->n_dofs_per_cell());
3476 AssertIndexRange(components.first_tensor_component, this->n_components());
3477 AssertIndexRange(components.first_tensor_component + dim,
3478 this->n_components() + 1);
3479
3480 if (is_primitive(shape_function))
3481 return ((system_to_component_table[shape_function].first >=
3482 components.first_tensor_component) &&
3483 (system_to_component_table[shape_function].first <
3484 components.first_tensor_component +
3486 else
3487 // Return whether there is overlap between the nonzero components
3488 // of the current shape function and the selected components of
3489 // the extractor:
3490 {
3491 for (unsigned int i = components.first_tensor_component;
3492 i < components.first_tensor_component +
3494 ++i)
3495 if (nonzero_components[shape_function][i])
3496 return true;
3497 return false;
3498 }
3499}
3500
3501
3502
3503template <int dim, int spacedim>
3504inline bool
3506 const unsigned int index) const
3507{
3508 AssertIndexRange(index, this->n_dofs_per_cell());
3509 return restriction_is_additive_flags[index];
3510}
3511
3512
3513
3514template <int dim, int spacedim>
3515inline const ComponentMask &
3517{
3519 return nonzero_components[i];
3520}
3521
3522
3523
3524template <int dim, int spacedim>
3525inline unsigned int
3527{
3529 return n_nonzero_components_table[i];
3530}
3531
3532
3533
3534template <int dim, int spacedim>
3535inline bool
3537{
3538 return cached_primitivity;
3539}
3540
3541
3542
3543template <int dim, int spacedim>
3544inline bool
3545FiniteElement<dim, spacedim>::is_primitive(const unsigned int i) const
3546{
3548
3549 // return primitivity of a shape
3550 // function by checking whether it
3551 // has more than one non-zero
3552 // component or not. we could cache
3553 // this value in an array of bools,
3554 // but accessing a bit-vector (as
3555 // std::vector<bool> is) is
3556 // probably more expensive than
3557 // just comparing against 1
3558 //
3559 // for good measure, short circuit the test
3560 // if the entire FE is primitive
3561 return (is_primitive() || (n_nonzero_components_table[i] == 1));
3562}
3563
3564
3565
3566template <int dim, int spacedim>
3567inline const Table<2, bool> &
3569{
3571}
3572
3573
3574
3575template <int dim, int spacedim>
3576inline GeometryPrimitive
3578 const unsigned int cell_dof_index) const
3579{
3580 AssertIndexRange(cell_dof_index, this->n_dofs_per_cell());
3581
3582 // just go through the usual cases, taking into account how DoFs
3583 // are enumerated on the reference cell
3584 if (cell_dof_index < this->get_first_line_index())
3586 else if (cell_dof_index < this->get_first_quad_index(0))
3588 else if (cell_dof_index < this->get_first_hex_index())
3590 else
3592}
3593
3594#endif
3595
3597
3598#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 types::geometric_orientation combined_orientation=numbers::default_geometric_orientation) const override
unsigned int get_first_line_index() const
const unsigned int components
Definition fe_data.h:446
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
bool shape_function_belongs_to(const unsigned int shape_function, const FEValuesExtractors::Vector &components) 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:2707
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:2568
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 shape_function_belongs_to(const unsigned int shape_function, const FEValuesExtractors::Scalar &component) const
bool has_face_support_points(const unsigned int face_no=0) const
const bool cached_primitivity
Definition fe.h:2714
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:2523
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:2597
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
bool shape_function_belongs_to(const unsigned int shape_function, const FEValuesExtractors::SymmetricTensor< 2 > &components) 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:2659
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:2580
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:2610
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 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
unsigned int adjust_quad_dof_index_for_face_orientation(const unsigned int index, const unsigned int face_no, const types::geometric_orientation combined_orientation) 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:2646
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
bool shape_function_belongs_to(const unsigned int shape_function, const FEValuesExtractors::Tensor< 2 > &components) 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:2615
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
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const types::geometric_orientation combined_orientation) const
const std::vector< bool > restriction_is_additive_flags
Definition fe.h:2689
virtual ~FiniteElement() override=default
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const types::geometric_orientation combined_orientation=numbers::default_geometric_orientation) const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition fe.h:2682
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:2653
std::vector< Point< dim > > unit_support_points
Definition fe.h:2561
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:2627
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:2549
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
Table< 2, bool > local_dof_sparsity_pattern
Definition fe.h:2722
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 const Table< 2, bool > & get_local_dof_sparsity_pattern() 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:2574
unsigned int component_to_system_index(const unsigned int component, const unsigned int index) const
const std::vector< ComponentMask > nonzero_components
Definition fe.h:2698
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:2537
types::global_dof_index first_block_of_base(const unsigned int b) const
Abstract base class for mapping classes.
Definition mapping.h:320
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
Point< 2 > first
Definition grid_out.cc:4632
#define DeclException0(Exception0)
static ::ExceptionBase & ExcShapeFunctionNotPrimitive(int arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
static ::ExceptionBase & ExcProjectionVoid()
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcEmbeddingVoid()
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcInterpolationNotImplemented()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcFEHasNoSupportPoints()
UpdateFlags
std::size_t size
Definition mpi.cc:745
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition utilities.cc:466
constexpr types::geometric_orientation default_geometric_orientation
Definition types.h:346
unsigned char geometric_orientation
Definition types.h:40