deal.II version GIT relicensing-3535-g02565b63bb 2025-06-19 13:50: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
dof_tools.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1999 - 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_dof_tools_h
16#define dealii_dof_tools_h
17
18
19#include <deal.II/base/config.h>
20
23#include <deal.II/base/point.h>
24
26
28
31
32#include <map>
33#include <ostream>
34#include <set>
35#include <vector>
36
37
39
40// Forward declarations
41#ifndef DOXYGEN
42class BlockMask;
43template <int dim, typename RangeNumberType>
44class Function;
45template <int dim, int spacedim>
46class FiniteElement;
47namespace hp
48{
49 template <int dim, int spacedim>
50 class MappingCollection;
51 template <int dim, int spacedim>
52 class FECollection;
53} // namespace hp
54template <typename MeshType>
56class InterGridMap;
57template <int dim, int spacedim>
58class Mapping;
59template <int dim, class T>
60class Table;
61template <typename Number>
62class Vector;
63
64namespace GridTools
65{
66 template <typename CellIterator>
67 struct PeriodicFacePair;
68}
69
70namespace DoFTools
71{
72 namespace internal
73 {
74 /*
75 * Default value of the face_has_flux_coupling parameter of
76 * make_flux_sparsity_pattern. Defined here (instead of using a default
77 * lambda in the parameter list) to avoid a bug in gcc where the same lambda
78 * gets defined multiple times.
79 */
80 template <int dim, int spacedim>
81 inline bool
82 always_couple_on_faces(
84 const unsigned int)
85 {
86 return true;
87 }
88 } // namespace internal
89} // namespace DoFTools
90
91#endif
92
225namespace DoFTools
226{
257
270 template <int dim, int spacedim>
271 void
273 const Table<2, Coupling> &table_by_component,
274 std::vector<Table<2, Coupling>> &tables_by_block);
275
281 template <int dim, int spacedim>
285 const Table<2, Coupling> &component_couplings);
286
294 template <int dim, int spacedim>
295 std::vector<Table<2, Coupling>>
298 const Table<2, Coupling> &component_couplings);
422 template <int dim, int spacedim, typename number = double>
423 void
425 const DoFHandler<dim, spacedim> &dof_handler,
426 SparsityPatternBase &sparsity_pattern,
427 const AffineConstraints<number> &constraints = {},
428 const bool keep_constrained_dofs = true,
430
496 template <int dim, int spacedim, typename number = double>
497 void
499 const DoFHandler<dim, spacedim> &dof_handler,
500 const Table<2, Coupling> &coupling,
501 SparsityPatternBase &sparsity_pattern,
502 const AffineConstraints<number> &constraints = {},
503 const bool keep_constrained_dofs = true,
505
526 template <int dim, int spacedim>
527 void
529 const DoFHandler<dim, spacedim> &dof_col,
530 SparsityPatternBase &sparsity);
531
577 template <int dim, int spacedim>
578 void
580 SparsityPatternBase &sparsity_pattern);
581
590 template <int dim, int spacedim, typename number>
591 void
593 const DoFHandler<dim, spacedim> &dof_handler,
594 SparsityPatternBase &sparsity_pattern,
595 const AffineConstraints<number> &constraints,
596 const bool keep_constrained_dofs = true,
598
599
619 template <int dim, int spacedim>
620 void
622 const DoFHandler<dim, spacedim> &dof,
623 SparsityPatternBase &sparsity,
624 const Table<2, Coupling> &cell_integrals_mask,
625 const Table<2, Coupling> &face_integrals_mask,
627
628
655 template <int dim, int spacedim, typename number>
656 void
658 const DoFHandler<dim, spacedim> &dof,
659 SparsityPatternBase &sparsity,
660 const AffineConstraints<number> &constraints,
661 const bool keep_constrained_dofs,
662 const Table<2, Coupling> &couplings,
663 const Table<2, Coupling> &face_couplings,
664 const types::subdomain_id subdomain_id,
665 const std::function<
667 const unsigned int)> &face_has_flux_coupling =
668 &internal::always_couple_on_faces<dim, spacedim>);
669
679 template <int dim, int spacedim>
680 void
682 const DoFHandler<dim, spacedim> &dof,
683 const std::vector<types::global_dof_index> &dof_to_boundary_mapping,
684 SparsityPatternBase &sparsity_pattern);
685
703 template <int dim, int spacedim, typename number>
704 void
706 const DoFHandler<dim, spacedim> &dof,
707 const std::map<types::boundary_id, const Function<spacedim, number> *>
708 &boundary_ids,
709 const std::vector<types::global_dof_index> &dof_to_boundary_mapping,
710 SparsityPatternBase &sparsity);
711
759 template <int dim, int spacedim, typename number>
760 void
762 AffineConstraints<number> &constraints);
763
831 template <int dim, int spacedim>
832 void
834 const DoFHandler<dim, spacedim> &coarse_grid,
835 const unsigned int coarse_component,
836 const DoFHandler<dim, spacedim> &fine_grid,
837 const unsigned int fine_component,
838 const InterGridMap<DoFHandler<dim, spacedim>> &coarse_to_fine_grid_map,
839 AffineConstraints<double> &constraints);
840
841
858 template <int dim, int spacedim>
859 void
861 const DoFHandler<dim, spacedim> &coarse_grid,
862 const unsigned int coarse_component,
863 const DoFHandler<dim, spacedim> &fine_grid,
864 const unsigned int fine_component,
865 const InterGridMap<DoFHandler<dim, spacedim>> &coarse_to_fine_grid_map,
866 std::vector<std::map<types::global_dof_index, float>>
867 &transfer_representation);
868
955 template <typename FaceIterator, typename number>
956 void
958 const FaceIterator &face_1,
960 AffineConstraints<number> &constraints,
961 const ComponentMask &component_mask = {},
962 const types::geometric_orientation combined_orientation =
965 const std::vector<unsigned int> &first_vector_components =
966 std::vector<unsigned int>(),
967 const number periodicity_factor = 1.);
968
969
989 template <int dim, int spacedim, typename number>
990 void
992 const std::vector<GridTools::PeriodicFacePair<
993 typename DoFHandler<dim, spacedim>::cell_iterator>> &periodic_faces,
994 AffineConstraints<number> &constraints,
995 const ComponentMask &component_mask = {},
996 const std::vector<unsigned int> &first_vector_components =
997 std::vector<unsigned int>(),
998 const number periodicity_factor = 1.);
999
1000
1001
1031 template <int dim, int spacedim, typename number>
1032 void
1034 const types::boundary_id b_id1,
1035 const types::boundary_id b_id2,
1036 const unsigned int direction,
1037 AffineConstraints<number> &constraints,
1038 const ComponentMask &component_mask = {},
1039 const number periodicity_factor = 1.);
1040
1041
1042
1064 template <int dim, int spacedim, typename number>
1065 void
1067 const types::boundary_id b_id,
1068 const unsigned int direction,
1069 AffineConstraints<number> &constraints,
1070 const ComponentMask &component_mask = {},
1071 const number periodicity_factor = 1.);
1072
1073 namespace internal
1074 {
1101 template <typename FaceIterator, typename number>
1102 void
1104 const FaceIterator &face_1,
1106 const FullMatrix<double> &transformation,
1107 AffineConstraints<number> &affine_constraints,
1108 const ComponentMask &component_mask,
1109 const types::geometric_orientation combined_orientation,
1110 const number periodicity_factor,
1111 const unsigned int level = numbers::invalid_unsigned_int);
1112 } // namespace internal
1113
1131 template <int dim, int spacedim>
1132 IndexSet
1134
1170 template <int dim, int spacedim>
1171 IndexSet
1172 extract_dofs(const DoFHandler<dim, spacedim> &dof_handler,
1173 const ComponentMask &component_mask);
1174
1198 template <int dim, int spacedim>
1199 IndexSet
1200 extract_dofs(const DoFHandler<dim, spacedim> &dof_handler,
1201 const BlockMask &block_mask);
1202
1207 template <int dim, int spacedim>
1208 void
1209 extract_level_dofs(const unsigned int level,
1210 const DoFHandler<dim, spacedim> &dof,
1211 const ComponentMask &component_mask,
1212 std::vector<bool> &selected_dofs);
1213
1218 template <int dim, int spacedim>
1219 void
1220 extract_level_dofs(const unsigned int level,
1221 const DoFHandler<dim, spacedim> &dof,
1222 const BlockMask &component_mask,
1223 std::vector<bool> &selected_dofs);
1224
1267 template <int dim, int spacedim>
1268 IndexSet
1270 const ComponentMask &component_mask = {},
1271 const std::set<types::boundary_id> &boundary_ids = {});
1272
1289 template <int dim, int spacedim>
1290 void
1292 const DoFHandler<dim, spacedim> &dof_handler,
1293 const ComponentMask &component_mask,
1294 std::vector<bool> &selected_dofs,
1295 const std::set<types::boundary_id> &boundary_ids =
1296 std::set<types::boundary_id>());
1297
1325 template <int dim, int spacedim, typename number = double>
1326 IndexSet
1328 const DoFHandler<dim, spacedim> &dof_handler,
1329 const std::function<
1331 &predicate,
1332 const AffineConstraints<number> &constraints = {});
1333
1369 template <int dim, int spacedim>
1370 std::vector<std::vector<bool>>
1372 const ComponentMask &component_mask = {});
1373
1379 template <int dim, int spacedim>
1380 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
1381 "Use the other function that returns the constant modes by value, rather than via an argument.")
1382 void extract_constant_modes(const DoFHandler<dim, spacedim> &dof_handler,
1383 const ComponentMask &component_mask,
1384 std::vector<std::vector<bool>> &constant_modes);
1385
1389 template <int dim, int spacedim>
1390 std::vector<std::vector<bool>>
1391 extract_level_constant_modes(const unsigned int level,
1392 const DoFHandler<dim, spacedim> &dof_handler,
1393 const ComponentMask &component_mask = {});
1394
1400 template <int dim, int spacedim>
1401 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
1402 "Use the other function that returns the constant modes by value, rather than via an argument.")
1404 const unsigned int level,
1405 const DoFHandler<dim, spacedim> &dof_handler,
1406 const ComponentMask &component_mask,
1407 std::vector<std::vector<bool>> &constant_modes);
1408
1414 template <int dim, int spacedim>
1415 std::vector<std::vector<double>>
1416 extract_rigid_body_modes(const Mapping<dim, spacedim> &mapping,
1417 const DoFHandler<dim, spacedim> &dof_handler,
1418 const ComponentMask &component_mask = {});
1419
1423 template <int dim, int spacedim>
1424 std::vector<std::vector<double>>
1425 extract_level_rigid_body_modes(const unsigned int level,
1426 const Mapping<dim, spacedim> &mapping,
1427 const DoFHandler<dim, spacedim> &dof_handler,
1428 const ComponentMask &component_mask = {});
1515 template <int dim, int spacedim>
1516 std::map<typename DoFHandler<dim - 1, spacedim>::active_cell_iterator,
1517 std::pair<typename DoFHandler<dim, spacedim>::active_cell_iterator,
1518 unsigned int>>
1522 &c1_to_c0,
1523 const DoFHandler<dim, spacedim> &c0_dh,
1524 const DoFHandler<dim - 1, spacedim> &c1_dh);
1525
1543 template <int dim, int spacedim>
1544 void
1546 const types::subdomain_id subdomain_id,
1547 std::vector<bool> &selected_dofs);
1548
1563 template <int dim, int spacedim>
1564 IndexSet
1566
1583 template <int dim, int spacedim>
1586 IndexSet &dof_set);
1587
1594 template <int dim, int spacedim>
1595 IndexSet
1597 const DoFHandler<dim, spacedim> &dof_handler,
1598 const unsigned int level);
1599
1608 template <int dim, int spacedim>
1611 const DoFHandler<dim, spacedim> &dof_handler,
1612 IndexSet &dof_set,
1613 const unsigned int level);
1614
1624 template <int dim, int spacedim>
1625 IndexSet
1627
1638 template <int dim, int spacedim>
1641 IndexSet &dof_set);
1642
1652 template <int dim, int spacedim>
1653 std::vector<IndexSet>
1655 const ComponentMask &components = {});
1656
1670 template <int dim, int spacedim>
1671 std::vector<IndexSet>
1673 const DoFHandler<dim, spacedim> &dof_handler);
1674
1689 template <int dim, int spacedim>
1690 std::vector<IndexSet>
1692 const DoFHandler<dim, spacedim> &dof_handler);
1693
1698 template <int dim, int spacedim>
1699 IndexSet
1701 const DoFHandler<dim, spacedim> &dof_handler,
1702 const unsigned int level);
1703
1710 template <int dim, int spacedim>
1713 const DoFHandler<dim, spacedim> &dof_handler,
1714 const unsigned int level,
1715 IndexSet &dof_set);
1716
1717
1748 template <int dim, int spacedim>
1749 void
1751 std::vector<types::subdomain_id> &subdomain);
1752
1778 template <int dim, int spacedim>
1779 unsigned int
1781 const DoFHandler<dim, spacedim> &dof_handler,
1782 const types::subdomain_id subdomain);
1783
1803 template <int dim, int spacedim>
1804 void
1806 const DoFHandler<dim, spacedim> &dof_handler,
1807 const types::subdomain_id subdomain,
1808 std::vector<unsigned int> &n_dofs_on_subdomain);
1809
1831 template <int dim, int spacedim>
1832 IndexSet
1834 const DoFHandler<dim, spacedim> &dof_handler,
1835 const types::subdomain_id subdomain);
1892 template <int dim, int spacedim>
1893 std::vector<types::global_dof_index>
1895 const std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>
1896 &patch);
1897
1918 template <int dim, int spacedim>
1919 void
1921 const DoFHandler<dim, spacedim> &dof_handler,
1922 const unsigned int level,
1923 const std::vector<bool> &selected_dofs = {},
1924 const types::global_dof_index offset = 0);
1925
1977 template <int dim, int spacedim>
1978 std::vector<unsigned int>
1980 const DoFHandler<dim, spacedim> &dof_handler,
1981 const unsigned int level,
1982 const bool interior_dofs_only,
1983 const bool boundary_patches = false,
1984 const bool level_boundary_patches = false,
1985 const bool single_cell_patches = false,
1986 const bool invert_vertex_mapping = false);
1987
2002 template <int dim, int spacedim>
2003 std::vector<unsigned int>
2005 const DoFHandler<dim, spacedim> &dof_handler,
2006 const unsigned int level,
2007 const BlockMask &exclude_boundary_dofs = BlockMask(),
2008 const bool boundary_patches = false,
2009 const bool level_boundary_patches = false,
2010 const bool single_cell_patches = false,
2011 const bool invert_vertex_mapping = false);
2012
2050 template <int dim, int spacedim>
2051 void
2053 const DoFHandler<dim, spacedim> &dof_handler,
2054 const unsigned int level,
2055 const bool interior_dofs_only,
2056 const bool boundary_dofs = false);
2057
2077 template <int dim, int spacedim>
2078 void
2080 const DoFHandler<dim, spacedim> &dof_handler,
2081 const unsigned int level,
2082 const bool interior_dofs_only = false);
2083
2126 template <int dim, int spacedim>
2127 std::vector<types::global_dof_index>
2129 const DoFHandler<dim, spacedim> &dof_handler,
2130 const bool vector_valued_once = false,
2131 const std::vector<unsigned int> &target_component = {});
2132
2149 template <int dim, int spacedim>
2150 std::vector<types::global_dof_index>
2152 const std::vector<unsigned int> &target_block =
2153 std::vector<unsigned int>());
2154
2170 template <int dim, int spacedim>
2173 std::vector<unsigned int> &active_fe_indices);
2174
2202 template <int dim, int spacedim>
2203 unsigned int
2205 const std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>
2206 &patch);
2207
2230 template <int dim, int spacedim>
2231 void
2233 std::vector<types::global_dof_index> &mapping);
2234
2245 template <int dim, int spacedim>
2246 void
2248 const std::set<types::boundary_id> &boundary_ids,
2249 std::vector<types::global_dof_index> &mapping);
2250
2282 template <int dim, int spacedim>
2283 void
2285 const DoFHandler<dim, spacedim> &dof_handler,
2286 std::vector<Point<spacedim>> &support_points,
2287 const ComponentMask &mask = {});
2288
2292 template <int dim, int spacedim>
2293 void
2296 const DoFHandler<dim, spacedim> &dof_handler,
2297 std::vector<Point<spacedim>> &support_points,
2298 const ComponentMask &mask = {});
2299
2329 template <int dim, int spacedim>
2330 std::map<types::global_dof_index, Point<spacedim>>
2332 const DoFHandler<dim, spacedim> &dof_handler,
2333 const ComponentMask &mask = {});
2334
2338 template <int dim, int spacedim>
2339 std::map<types::global_dof_index, Point<spacedim>>
2342 const DoFHandler<dim, spacedim> &dof_handler,
2343 const ComponentMask &mask = {});
2344
2350 template <int dim, int spacedim>
2353 const Mapping<dim, spacedim> &mapping,
2354 const DoFHandler<dim, spacedim> &dof_handler,
2355 std::map<types::global_dof_index, Point<spacedim>> &support_points,
2356 const ComponentMask &mask = {});
2357
2363 template <int dim, int spacedim>
2367 const DoFHandler<dim, spacedim> &dof_handler,
2368 std::map<types::global_dof_index, Point<spacedim>> &support_points,
2369 const ComponentMask &mask = {});
2370
2371
2392 template <int dim, int spacedim, class Comp>
2393 void
2395 const Mapping<dim, spacedim> &mapping,
2396 const DoFHandler<dim, spacedim> &dof_handler,
2398 &point_to_index_map);
2436 template <int dim, int spacedim, typename Number>
2437 void
2439 const Vector<Number> &cell_data,
2440 Vector<double> &dof_data,
2441 const unsigned int component = 0);
2442
2443
2521 template <int spacedim>
2522 void
2524 std::ostream &out,
2525 const std::map<types::global_dof_index, Point<spacedim>> &support_points);
2526
2527
2568 template <int dim, int spacedim, typename number>
2569 void
2571 const DoFHandler<dim, spacedim> &dof,
2572 const types::boundary_id boundary_id,
2573 AffineConstraints<number> &zero_boundary_constraints,
2574 const ComponentMask &component_mask = {});
2575
2586 template <int dim, int spacedim, typename number>
2587 void
2589 const DoFHandler<dim, spacedim> &dof,
2590 AffineConstraints<number> &zero_boundary_constraints,
2591 const ComponentMask &component_mask = {});
2592
2637} // namespace DoFTools
2638
2639
2640
2641/* ------------------------- inline functions -------------- */
2642
2643#ifndef DOXYGEN
2644
2645namespace DoFTools
2646{
2652 inline Coupling
2653 operator|=(Coupling &c1, const Coupling c2)
2654 {
2655 if (c2 == always)
2656 c1 = always;
2657 else if (c1 != always && c2 == nonzero)
2658 return c1 = nonzero;
2659 return c1;
2660 }
2661
2662
2668 inline Coupling
2669 operator|(const Coupling c1, const Coupling c2)
2670 {
2671 if (c1 == always || c2 == always)
2672 return always;
2673 if (c1 == nonzero || c2 == nonzero)
2674 return nonzero;
2675 return none;
2676 }
2677
2678
2679 // ---------------------- inline and template functions --------------------
2680
2681 template <int dim, int spacedim, class Comp>
2682 void
2684 const Mapping<dim, spacedim> &mapping,
2685 const DoFHandler<dim, spacedim> &dof_handler,
2687 &point_to_index_map)
2688 {
2689 // let the checking of arguments be
2690 // done by the function first
2691 // called
2692 std::vector<Point<spacedim>> support_points(dof_handler.n_dofs());
2693 map_dofs_to_support_points(mapping, dof_handler, support_points);
2694 // now copy over the results of the
2695 // previous function into the
2696 // output arg
2697 point_to_index_map.clear();
2698 for (types::global_dof_index i = 0; i < dof_handler.n_dofs(); ++i)
2699 point_to_index_map[support_points[i]] = i;
2700 }
2701} // namespace DoFTools
2702
2703#endif
2704
2706
2707#endif
types::global_dof_index n_dofs() const
Abstract base class for mapping classes.
Definition mapping.h:320
Definition point.h:113
#define DEAL_II_DEPRECATED
Definition config.h:281
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:243
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
unsigned int level
Definition grid_out.cc:4635
#define DeclException0(Exception0)
static ::ExceptionBase & ExcGridNotCoarser()
static ::ExceptionBase & ExcInvalidBoundaryIndicator()
static ::ExceptionBase & ExcFiniteElementsDontMatch()
static ::ExceptionBase & ExcNoFESelected()
static ::ExceptionBase & ExcGridsDontMatch()
typename ActiveSelector::cell_iterator cell_iterator
typename ActiveSelector::active_cell_iterator active_cell_iterator
void compute_intergrid_transfer_representation(const DoFHandler< dim, spacedim > &coarse_grid, const unsigned int coarse_component, const DoFHandler< dim, spacedim > &fine_grid, const unsigned int fine_component, const InterGridMap< DoFHandler< dim, spacedim > > &coarse_to_fine_grid_map, std::vector< std::map< types::global_dof_index, float > > &transfer_representation)
void make_hanging_node_constraints(const DoFHandler< dim, spacedim > &dof_handler, AffineConstraints< number > &constraints)
void make_sparsity_pattern(const DoFHandler< dim, spacedim > &dof_handler, SparsityPatternBase &sparsity_pattern, const AffineConstraints< number > &constraints={}, const bool keep_constrained_dofs=true, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id)
void compute_intergrid_constraints(const DoFHandler< dim, spacedim > &coarse_grid, const unsigned int coarse_component, const DoFHandler< dim, spacedim > &fine_grid, const unsigned int fine_component, const InterGridMap< DoFHandler< dim, spacedim > > &coarse_to_fine_grid_map, AffineConstraints< double > &constraints)
void make_zero_boundary_constraints(const DoFHandler< dim, spacedim > &dof, const types::boundary_id boundary_id, AffineConstraints< number > &zero_boundary_constraints, const ComponentMask &component_mask={})
void make_flux_sparsity_pattern(const DoFHandler< dim, spacedim > &dof_handler, SparsityPatternBase &sparsity_pattern)
void make_boundary_sparsity_pattern(const DoFHandler< dim, spacedim > &dof, const std::vector< types::global_dof_index > &dof_to_boundary_mapping, SparsityPatternBase &sparsity_pattern)
UpdateFlags & operator|=(UpdateFlags &f1, const UpdateFlags f2)
InitializeLibrary operator|(const InitializeLibrary f1, const InitializeLibrary f2)
void set_periodicity_constraints(const FaceIterator &face_1, const std_cxx20::type_identity_t< FaceIterator > &face_2, const FullMatrix< double > &transformation, AffineConstraints< number > &affine_constraints, const ComponentMask &component_mask, const types::geometric_orientation combined_orientation, const number periodicity_factor, const unsigned int level=numbers::invalid_unsigned_int)
void get_subdomain_association(const DoFHandler< dim, spacedim > &dof_handler, std::vector< types::subdomain_id > &subdomain)
IndexSet dof_indices_with_subdomain_association(const DoFHandler< dim, spacedim > &dof_handler, const types::subdomain_id subdomain)
IndexSet extract_dofs_with_support_contained_within(const DoFHandler< dim, spacedim > &dof_handler, const std::function< bool(const typename DoFHandler< dim, spacedim >::active_cell_iterator &)> &predicate, const AffineConstraints< number > &constraints={})
Definition dof_tools.cc:826
std::vector< IndexSet > locally_owned_dofs_per_subdomain(const DoFHandler< dim, spacedim > &dof_handler)
IndexSet extract_boundary_dofs(const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &component_mask={}, const std::set< types::boundary_id > &boundary_ids={})
Definition dof_tools.cc:616
std::vector< std::vector< bool > > extract_level_constant_modes(const unsigned int level, const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &component_mask={})
IndexSet extract_locally_relevant_dofs(const DoFHandler< dim, spacedim > &dof_handler)
std::vector< std::vector< bool > > extract_constant_modes(const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &component_mask={})
void map_support_points_to_dofs(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof_handler, std::map< Point< spacedim >, types::global_dof_index, Comp > &point_to_index_map)
std::map< typename DoFHandler< dim - 1, spacedim >::active_cell_iterator, std::pair< typename DoFHandler< dim, spacedim >::active_cell_iterator, unsigned int > > map_boundary_to_bulk_dof_iterators(const std::map< typename Triangulation< dim - 1, spacedim >::cell_iterator, typename Triangulation< dim, spacedim >::face_iterator > &c1_to_c0, const DoFHandler< dim, spacedim > &c0_dh, const DoFHandler< dim - 1, spacedim > &c1_dh)
void get_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, std::vector< unsigned int > &active_fe_indices)
void make_cell_patches(SparsityPattern &block_list, const DoFHandler< dim, spacedim > &dof_handler, const unsigned int level, const std::vector< bool > &selected_dofs={}, const types::global_dof_index offset=0)
IndexSet extract_dofs(const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &component_mask)
Definition dof_tools.cc:416
IndexSet extract_locally_active_dofs(const DoFHandler< dim, spacedim > &dof_handler)
IndexSet extract_locally_relevant_level_dofs(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int level)
std::vector< std::vector< double > > extract_rigid_body_modes(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &component_mask={})
void map_dofs_to_support_points(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof_handler, std::vector< Point< spacedim > > &support_points, const ComponentMask &mask={})
void write_gnuplot_dof_support_point_info(std::ostream &out, const std::map< types::global_dof_index, Point< spacedim > > &support_points)
void extract_subdomain_dofs(const DoFHandler< dim, spacedim > &dof_handler, const types::subdomain_id subdomain_id, std::vector< bool > &selected_dofs)
IndexSet extract_locally_active_level_dofs(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int level)
std::vector< unsigned int > make_vertex_patches(SparsityPattern &block_list, const DoFHandler< dim, spacedim > &dof_handler, const unsigned int level, const bool interior_dofs_only, const bool boundary_patches=false, const bool level_boundary_patches=false, const bool single_cell_patches=false, const bool invert_vertex_mapping=false)
std::vector< types::global_dof_index > count_dofs_per_fe_block(const DoFHandler< dim, spacedim > &dof, const std::vector< unsigned int > &target_block=std::vector< unsigned int >())
void extract_level_dofs(const unsigned int level, const DoFHandler< dim, spacedim > &dof, const ComponentMask &component_mask, std::vector< bool > &selected_dofs)
Definition dof_tools.cc:503
void extract_dofs_with_support_on_boundary(const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &component_mask, std::vector< bool > &selected_dofs, const std::set< types::boundary_id > &boundary_ids=std::set< types::boundary_id >())
Definition dof_tools.cc:744
void distribute_cell_to_dof_vector(const DoFHandler< dim, spacedim > &dof_handler, const Vector< Number > &cell_data, Vector< double > &dof_data, const unsigned int component=0)
Definition dof_tools.cc:332
std::vector< types::global_dof_index > count_dofs_per_fe_component(const DoFHandler< dim, spacedim > &dof_handler, const bool vector_valued_once=false, const std::vector< unsigned int > &target_component={})
void make_child_patches(SparsityPattern &block_list, const DoFHandler< dim, spacedim > &dof_handler, const unsigned int level, const bool interior_dofs_only, const bool boundary_dofs=false)
void map_dof_to_boundary_indices(const DoFHandler< dim, spacedim > &dof_handler, std::vector< types::global_dof_index > &mapping)
std::vector< IndexSet > locally_relevant_dofs_per_subdomain(const DoFHandler< dim, spacedim > &dof_handler)
void make_single_patch(SparsityPattern &block_list, const DoFHandler< dim, spacedim > &dof_handler, const unsigned int level, const bool interior_dofs_only=false)
std::vector< types::global_dof_index > get_dofs_on_patch(const std::vector< typename DoFHandler< dim, spacedim >::active_cell_iterator > &patch)
std::vector< IndexSet > locally_owned_dofs_per_component(const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &components={})
Definition dof_tools.cc:469
unsigned int count_dofs_with_subdomain_association(const DoFHandler< dim, spacedim > &dof_handler, const types::subdomain_id subdomain)
unsigned int count_dofs_on_patch(const std::vector< typename DoFHandler< dim, spacedim >::active_cell_iterator > &patch)
void convert_couplings_to_blocks(const DoFHandler< dim, spacedim > &dof_handler, const Table< 2, Coupling > &table_by_component, std::vector< Table< 2, Coupling > > &tables_by_block)
Table< 2, Coupling > dof_couplings_from_component_couplings(const FiniteElement< dim, spacedim > &fe, const Table< 2, Coupling > &component_couplings)
IndexSet extract_hanging_node_dofs(const DoFHandler< dim, spacedim > &dof_handler)
void make_periodicity_constraints(const FaceIterator &face_1, const std_cxx20::type_identity_t< FaceIterator > &face_2, AffineConstraints< number > &constraints, const ComponentMask &component_mask={}, const types::geometric_orientation combined_orientation=numbers::default_geometric_orientation, const FullMatrix< double > &matrix=FullMatrix< double >(), const std::vector< unsigned int > &first_vector_components=std::vector< unsigned int >(), const number periodicity_factor=1.)
std::vector< std::vector< double > > extract_level_rigid_body_modes(const unsigned int level, const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof_handler, const ComponentMask &component_mask={})
@ matrix
Contents is actually a matrix.
Definition hp.h:117
constexpr unsigned int invalid_unsigned_int
Definition types.h:238
constexpr types::subdomain_id invalid_subdomain_id
Definition types.h:381
constexpr types::geometric_orientation default_geometric_orientation
Definition types.h:352
typename type_identity< T >::type type_identity_t
Definition type_traits.h:95
STL namespace.
unsigned int subdomain_id
Definition types.h:52
unsigned char geometric_orientation
Definition types.h:40