Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-2684-gc61376a70f 2025-02-22 15:30: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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
fe_tools.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2000 - 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_tools_H
16#define dealii_fe_tools_H
17
18
19
20#include <deal.II/base/config.h>
21
26#include <deal.II/base/tensor.h>
27
29
31
33
34#include <memory>
35#include <string>
36#include <vector>
37
38
40
41// Forward declarations
42#ifndef DOXYGEN
43template <typename number>
44class FullMatrix;
45template <int dim>
46class Quadrature;
47template <int dim, int spacedim>
48class FiniteElement;
49template <int dim, int spacedim>
51class DoFHandler;
52template <int dim>
54template <typename number>
56#endif
57
58
77namespace FETools
78{
87 template <int dim, int spacedim = dim>
89 {
90 public:
94 virtual std::unique_ptr<FiniteElement<dim, spacedim>>
95 get(const unsigned int degree) const = 0;
96
102 virtual std::unique_ptr<FiniteElement<dim, spacedim>>
103 get(const Quadrature<1> &quad) const = 0;
104
108 virtual ~FEFactoryBase() override = default;
109 };
110
120 template <class FE>
121 class FEFactory : public FEFactoryBase<FE::dimension, FE::space_dimension>
122 {
123 public:
127 virtual std::unique_ptr<FiniteElement<FE::dimension, FE::space_dimension>>
128 get(const unsigned int degree) const override;
129
134 virtual std::unique_ptr<FiniteElement<FE::dimension, FE::space_dimension>>
135 get(const Quadrature<1> &quad) const override;
136 };
137
153 template <int dim, int spacedim>
154 void
156 std::vector<unsigned int> &renumbering,
157 std::vector<std::vector<unsigned int>> &start_indices);
158
174 template <int dim, int spacedim>
175 void
177 std::vector<types::global_dof_index> &renumbering,
178 std::vector<types::global_dof_index> &block_data,
179 bool return_start_indices = true);
180
194 template <int dim, typename number, int spacedim>
195 void
198 FullMatrix<number> &interpolation_matrix);
199
211 template <int dim, typename number, int spacedim>
212 void
215 FullMatrix<number> &interpolation_matrix);
216
228 template <int dim, typename number, int spacedim>
229 void
232 FullMatrix<number> &difference_matrix);
233
237 template <int dim, typename number, int spacedim>
238 void
241 FullMatrix<number> &matrix);
242
311 template <int dim, int spacedim>
314
353 template <int dim, typename number, int spacedim>
354 void
357 std::vector<std::vector<FullMatrix<number>>> &matrices,
358 const bool isotropic_only = false,
359 const double threshold = 1.e-12);
360
382 template <int dim, typename number, int spacedim>
383 void
385 const ArrayView<FullMatrix<number>> &matrices,
386 const unsigned int face_coarse,
387 const unsigned int face_fine,
388 const double threshold = 1.e-12);
389
421 template <int dim, typename number, int spacedim>
422 void
425 std::vector<std::vector<FullMatrix<number>>> &matrices,
426 const bool isotropic_only = false);
427
513 template <int dim, int spacedim>
514 void
517 const Quadrature<dim> &lhs_quadrature,
518 const Quadrature<dim> &rhs_quadrature,
520
528 template <int dim, int spacedim>
529 void
532 const Quadrature<dim> &quadrature,
533 FullMatrix<double> &I_q);
534
548 template <int dim>
549 void
551 const FullMatrix<double> &projection_matrix,
552 const std::vector<Tensor<1, dim>> &vector_of_tensors_at_qp,
553 std::vector<Tensor<1, dim>> &vector_of_tensors_at_nodes);
554
555
556
560 template <int dim>
561 void
563 const FullMatrix<double> &projection_matrix,
564 const std::vector<SymmetricTensor<2, dim>> &vector_of_tensors_at_qp,
565 std::vector<SymmetricTensor<2, dim>> &vector_of_tensors_at_nodes);
566
567
568
578 template <int dim, int spacedim>
579 void
582 const Quadrature<dim - 1> &lhs_quadrature,
583 const Quadrature<dim - 1> &rhs_quadrature,
585 const unsigned int face,
587
620 template <int dim, int spacedim = dim>
623
639 template <int dim, int spacedim, typename number>
640 void
642 const FiniteElement<dim, spacedim> &finite_element,
643 const std::vector<Vector<number>> &support_point_values,
644 std::vector<number> &dof_values);
645
646
647
685 template <int dim, int spacedim, class InVector, class OutVector>
686 void
688 const InVector &u1,
689 const DoFHandler<dim, spacedim> &dof2,
690 OutVector &u2);
691
711 template <int dim, int spacedim, class InVector, class OutVector>
712 void
714 const DoFHandler<dim, spacedim> &dof1,
715 const InVector &u1,
716 const DoFHandler<dim, spacedim> &dof2,
718 OutVector &u2);
719
733 template <int dim, class InVector, class OutVector, int spacedim>
734 void
736 const InVector &u1,
738 OutVector &u1_interpolated);
739
752 template <int dim, class InVector, class OutVector, int spacedim>
753 void
755 const DoFHandler<dim, spacedim> &dof1,
757 const InVector &u1,
758 const DoFHandler<dim, spacedim> &dof2,
760 OutVector &u1_interpolated);
761
771 template <int dim, class InVector, class OutVector, int spacedim>
772 void
774 const InVector &z1,
776 OutVector &z1_difference);
777
790 template <int dim, class InVector, class OutVector, int spacedim>
791 void
793 const DoFHandler<dim, spacedim> &dof1,
795 const InVector &z1,
796 const DoFHandler<dim, spacedim> &dof2,
798 OutVector &z1_difference);
799
800
801
810 template <int dim, class InVector, class OutVector, int spacedim>
811 void
813 const InVector &u1,
814 const DoFHandler<dim, spacedim> &dof2,
815 OutVector &u2);
816
877 template <int dim, class InVector, class OutVector, int spacedim>
878 void
880 const InVector &z1,
881 const DoFHandler<dim, spacedim> &dof2,
882 OutVector &z2);
883
896 template <int dim, class InVector, class OutVector, int spacedim>
897 void
899 const DoFHandler<dim, spacedim> &dof1,
900 const InVector &z1,
901 const DoFHandler<dim, spacedim> &dof2,
903 OutVector &z2);
904
920 template <int dim>
921 std::vector<unsigned int>
923
930 template <int dim>
931 std::vector<unsigned int>
933
1009 namespace Compositing
1010 {
1027 template <int dim, int spacedim>
1030 const std::vector<const FiniteElement<dim, spacedim> *> &fes,
1031 const std::vector<unsigned int> &multiplicities,
1032 const bool do_tensor_product = true);
1033
1039 template <int dim, int spacedim>
1042 const std::initializer_list<
1043 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>>
1044 &fe_systems);
1045
1052 template <int dim, int spacedim>
1055 const unsigned int N1,
1056 const FiniteElement<dim, spacedim> *fe2 = nullptr,
1057 const unsigned int N2 = 0,
1058 const FiniteElement<dim, spacedim> *fe3 = nullptr,
1059 const unsigned int N3 = 0,
1060 const FiniteElement<dim, spacedim> *fe4 = nullptr,
1061 const unsigned int N4 = 0,
1062 const FiniteElement<dim, spacedim> *fe5 = nullptr,
1063 const unsigned int N5 = 0);
1064
1077 template <int dim, int spacedim>
1078 std::vector<bool>
1080 const std::vector<const FiniteElement<dim, spacedim> *> &fes,
1081 const std::vector<unsigned int> &multiplicities);
1082
1088 template <int dim, int spacedim>
1089 std::vector<bool>
1091 const std::initializer_list<
1092 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>>
1093 &fe_systems);
1094
1112 template <int dim, int spacedim>
1113 DEAL_II_DEPRECATED std::vector<bool>
1116 const unsigned int N1,
1117 const FiniteElement<dim, spacedim> *fe2 = nullptr,
1118 const unsigned int N2 = 0,
1119 const FiniteElement<dim, spacedim> *fe3 = nullptr,
1120 const unsigned int N3 = 0,
1121 const FiniteElement<dim, spacedim> *fe4 = nullptr,
1122 const unsigned int N4 = 0,
1123 const FiniteElement<dim, spacedim> *fe5 = nullptr,
1124 const unsigned int N5 = 0);
1125
1126
1143 template <int dim, int spacedim>
1144 std::vector<ComponentMask>
1146 const std::vector<const FiniteElement<dim, spacedim> *> &fes,
1147 const std::vector<unsigned int> &multiplicities,
1148 const bool do_tensor_product = true);
1149
1155 template <int dim, int spacedim>
1156 std::vector<ComponentMask>
1158 const std::initializer_list<
1159 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>>
1160 &fe_systems);
1161
1183 template <int dim, int spacedim>
1184 DEAL_II_DEPRECATED std::vector<ComponentMask>
1187 const unsigned int N1,
1188 const FiniteElement<dim, spacedim> *fe2 = nullptr,
1189 const unsigned int N2 = 0,
1190 const FiniteElement<dim, spacedim> *fe3 = nullptr,
1191 const unsigned int N3 = 0,
1192 const FiniteElement<dim, spacedim> *fe4 = nullptr,
1193 const unsigned int N4 = 0,
1194 const FiniteElement<dim, spacedim> *fe5 = nullptr,
1195 const unsigned int N5 = 0,
1196 const bool do_tensor_product = true);
1197
1214 template <int dim, int spacedim>
1215 void
1217 std::vector<std::pair<std::pair<unsigned int, unsigned int>,
1218 unsigned int>> &system_to_base_table,
1219 std::vector<std::pair<unsigned int, unsigned int>>
1220 &system_to_component_table,
1221 std::vector<std::pair<std::pair<unsigned int, unsigned int>,
1222 unsigned int>> &component_to_base_table,
1223 const FiniteElement<dim, spacedim> &finite_element,
1224 const bool do_tensor_product = true);
1225
1241 template <int dim, int spacedim>
1242 void
1244 std::vector<std::pair<std::pair<unsigned int, unsigned int>,
1245 unsigned int>> &face_system_to_base_table,
1246 std::vector<std::pair<unsigned int, unsigned int>>
1247 &face_system_to_component_table,
1248 const FiniteElement<dim, spacedim> &finite_element,
1249 const bool do_tensor_product = true,
1250 const unsigned int face_no = 0 /*TODO*/);
1251
1252 } // namespace Compositing
1253
1254
1288 template <int dim, int spacedim = dim>
1289 std::unique_ptr<FiniteElement<dim, spacedim>>
1290 get_fe_by_name(const std::string &name);
1291
1334 template <int dim, int spacedim>
1335 void
1336 add_fe_name(const std::string &name,
1337 const FEFactoryBase<dim, spacedim> *factory);
1338
1349 std::string,
1350 << "Can't re-generate a finite element from the string '"
1351 << arg1 << "'.");
1352
1365 char,
1366 int,
1367 << "The dimension " << arg1
1368 << " in the finite element string must match "
1369 << "the space dimension " << arg2 << '.');
1370
1377
1391
1399 "You are using continuous elements on a grid with "
1400 "hanging nodes but without providing hanging node "
1401 "constraints. Use the respective function with "
1402 "additional AffineConstraints argument(s), instead.");
1415 int,
1416 int,
1417 int,
1418 int,
1419 << "This is a " << arg1 << 'x' << arg2 << " matrix, "
1420 << "but should be a " << arg3 << 'x' << arg4 << " matrix.");
1421
1428 double,
1429 << "Least squares fit leaves a gap of " << arg1);
1430
1437 int,
1438 int,
1439 << arg1 << " must be greater than " << arg2);
1440} // namespace FETools
1441
1442
1443#ifndef DOXYGEN
1444
1445namespace FETools
1446{
1447 template <class FE>
1448 std::unique_ptr<FiniteElement<FE::dimension, FE::space_dimension>>
1449 FEFactory<FE>::get(const unsigned int degree) const
1450 {
1451 return std::make_unique<FE>(degree);
1452 }
1453
1454 namespace Compositing
1455 {
1456 template <int dim, int spacedim>
1457 std::vector<bool>
1459 const std::initializer_list<
1460 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>>
1461 &fe_systems)
1462 {
1463 std::vector<const FiniteElement<dim, spacedim> *> fes;
1464 std::vector<unsigned int> multiplicities;
1465
1466 const auto extract =
1467 [&fes, &multiplicities](
1468 const std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
1469 unsigned int> &fe_system) {
1470 fes.push_back(fe_system.first.get());
1471 multiplicities.push_back(fe_system.second);
1472 };
1473
1474 for (const auto &p : fe_systems)
1475 extract(p);
1476
1477 return compute_restriction_is_additive_flags(fes, multiplicities);
1478 }
1479
1480
1481
1482 template <int dim, int spacedim>
1485 const std::initializer_list<
1486 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>>
1487 &fe_systems)
1488 {
1489 std::vector<const FiniteElement<dim, spacedim> *> fes;
1490 std::vector<unsigned int> multiplicities;
1491
1492 const auto extract =
1493 [&fes, &multiplicities](
1494 const std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
1495 unsigned int> &fe_system) {
1496 fes.push_back(fe_system.first.get());
1497 multiplicities.push_back(fe_system.second);
1498 };
1499
1500 for (const auto &p : fe_systems)
1501 extract(p);
1502
1503 return multiply_dof_numbers(fes, multiplicities, true);
1504 }
1505
1506
1507
1508 template <int dim, int spacedim>
1509 std::vector<ComponentMask>
1511 const std::initializer_list<
1512 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>>
1513 &fe_systems)
1514 {
1515 std::vector<const FiniteElement<dim, spacedim> *> fes;
1516 std::vector<unsigned int> multiplicities;
1517
1518 const auto extract =
1519 [&fes, &multiplicities](
1520 const std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
1521 unsigned int> &fe_system) {
1522 fes.push_back(fe_system.first.get());
1523 multiplicities.push_back(fe_system.second);
1524 };
1525
1526 for (const auto &p : fe_systems)
1527 extract(p);
1528
1529 return compute_nonzero_components(fes, multiplicities, true);
1530 }
1531 } // namespace Compositing
1532} // namespace FETools
1533
1534#endif
1535
1539
1540#endif /* dealii_fe_tools_H */
virtual ~FEFactoryBase() override=default
virtual std::unique_ptr< FiniteElement< dim, spacedim > > get(const Quadrature< 1 > &quad) const =0
virtual std::unique_ptr< FiniteElement< dim, spacedim > > get(const unsigned int degree) const =0
virtual std::unique_ptr< FiniteElement< FE::dimension, FE::space_dimension > > get(const Quadrature< 1 > &quad) const override
virtual std::unique_ptr< FiniteElement< FE::dimension, FE::space_dimension > > get(const unsigned int degree) const override
#define DEAL_II_DEPRECATED
Definition config.h:232
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:522
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:194
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:523
#define DeclException0(Exception0)
static ::ExceptionBase & ExcTriangulationMismatch()
static ::ExceptionBase & ExcInvalidFEName(std::string arg1)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
static ::ExceptionBase & ExcInvalidFEDimension(char arg1, int arg2)
static ::ExceptionBase & ExcLeastSquaresError(double arg1)
#define DeclException2(Exception2, type1, type2, outsequence)
static ::ExceptionBase & ExcFENotPrimitive()
static ::ExceptionBase & ExcHangingNodesNotAllowed()
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInvalidFE()
static ::ExceptionBase & ExcGridNotRefinedAtLeastOnce()
static ::ExceptionBase & ExcMatrixDimensionMismatch(int arg1, int arg2, int arg3, int arg4)
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcNotGreaterThan(int arg1, int arg2)
typename ActiveSelector::active_cell_iterator active_cell_iterator
void build_face_tables(std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > &face_system_to_base_table, std::vector< std::pair< unsigned int, unsigned int > > &face_system_to_component_table, const FiniteElement< dim, spacedim > &finite_element, const bool do_tensor_product=true, const unsigned int face_no=0)
std::vector< ComponentMask > compute_nonzero_components(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities, const bool do_tensor_product=true)
std::vector< bool > compute_restriction_is_additive_flags(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities)
void build_cell_tables(std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > &system_to_base_table, std::vector< std::pair< unsigned int, unsigned int > > &system_to_component_table, std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > &component_to_base_table, const FiniteElement< dim, spacedim > &finite_element, const bool do_tensor_product=true)
FiniteElementData< dim > multiply_dof_numbers(const std::vector< const FiniteElement< dim, spacedim > * > &fes, const std::vector< unsigned int > &multiplicities, const bool do_tensor_product=true)
void interpolation_difference(const DoFHandler< dim, spacedim > &dof1, const InVector &z1, const FiniteElement< dim, spacedim > &fe2, OutVector &z1_difference)
void compute_projection_from_quadrature_points(const FullMatrix< double > &projection_matrix, const std::vector< Tensor< 1, dim > > &vector_of_tensors_at_qp, std::vector< Tensor< 1, dim > > &vector_of_tensors_at_nodes)
void project_dg(const DoFHandler< dim, spacedim > &dof1, const InVector &u1, const DoFHandler< dim, spacedim > &dof2, OutVector &u2)
void get_back_interpolation_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &interpolation_matrix)
void get_projection_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &matrix)
void compute_component_wise(const FiniteElement< dim, spacedim > &fe, std::vector< unsigned int > &renumbering, std::vector< std::vector< unsigned int > > &start_indices)
void add_fe_name(const std::string &name, const FEFactoryBase< dim, spacedim > *factory)
std::vector< unsigned int > hierarchic_to_lexicographic_numbering(unsigned int degree)
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)
void compute_interpolation_to_quadrature_points_matrix(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, FullMatrix< double > &I_q)
Quadrature< dim > compute_nodal_quadrature(const FiniteElement< dim, spacedim > &fe)
Compute the nodal quadrature rule associated with an element.
std::unique_ptr< FiniteElement< dim, spacedim > > get_fe_by_name(const std::string &name)
void compute_projection_from_quadrature_points_matrix(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &lhs_quadrature, const Quadrature< dim > &rhs_quadrature, FullMatrix< double > &X)
void get_interpolation_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &interpolation_matrix)
void extrapolate(const DoFHandler< dim, spacedim > &dof1, const InVector &z1, const DoFHandler< dim, spacedim > &dof2, OutVector &z2)
void compute_embedding_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number > > > &matrices, const bool isotropic_only=false, const double threshold=1.e-12)
void compute_face_embedding_matrices(const FiniteElement< dim, spacedim > &fe, const ArrayView< FullMatrix< number > > &matrices, const unsigned int face_coarse, const unsigned int face_fine, const double threshold=1.e-12)
void interpolate(const DoFHandler< dim, spacedim > &dof1, const InVector &u1, const DoFHandler< dim, spacedim > &dof2, OutVector &u2)
void compute_block_renumbering(const FiniteElement< dim, spacedim > &fe, std::vector< types::global_dof_index > &renumbering, std::vector< types::global_dof_index > &block_data, bool return_start_indices=true)
std::vector< unsigned int > lexicographic_to_hierarchic_numbering(unsigned int degree)
void compute_projection_from_face_quadrature_points_matrix(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim - 1 > &lhs_quadrature, const Quadrature< dim - 1 > &rhs_quadrature, const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell, const unsigned int face, FullMatrix< double > &X)
FullMatrix< double > compute_node_matrix(const FiniteElement< dim, spacedim > &fe)
void get_interpolation_difference_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &difference_matrix)
void back_interpolate(const DoFHandler< dim, spacedim > &dof1, const InVector &u1, const FiniteElement< dim, spacedim > &fe2, OutVector &u1_interpolated)
void compute_projection_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number > > > &matrices, const bool isotropic_only=false)
constexpr ReturnType< rank, T >::value_type & extract(T &t, const ArrayType &indices)