Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3057-gb5b616642a 2025-04-10 20: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
dof_handler.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_dof_handler_h
16#define dealii_dof_handler_h
17
18
19
20#include <deal.II/base/config.h>
21
27#include <deal.II/base/types.h>
28
35
37
38#include <boost/serialization/split_member.hpp>
39#include <boost/signals2/connection.hpp>
40
41#include <map>
42#include <memory>
43#include <set>
44#include <vector>
45
47
48// Forward declarations
49#ifndef DOXYGEN
50template <int dim, int spacedim>
51class FiniteElement;
52template <int dim, int spacedim>
54class Triangulation;
55
56namespace internal
57{
58 namespace DoFHandlerImplementation
59 {
60 struct Implementation;
61
62 namespace Policy
63 {
64 template <int dim, int spacedim>
65 class PolicyBase;
66 struct Implementation;
67 } // namespace Policy
68 } // namespace DoFHandlerImplementation
69
70 namespace DoFAccessorImplementation
71 {
72 struct Implementation;
73 }
74
75 namespace DoFCellAccessorImplementation
76 {
77 struct Implementation;
78 }
79
80 namespace hp
81 {
82 namespace DoFHandlerImplementation
83 {
84 struct Implementation;
85 }
86 } // namespace hp
87} // namespace internal
88
89namespace parallel
90{
91 namespace distributed
92 {
93 template <int dim, int spacedim, typename VectorType>
94 class CellDataTransfer;
95 }
96} // namespace parallel
97
98template <int structdim, int dim, int spacedim, bool level_dof_access>
99class DoFAccessor;
100
101template <int dimension_, int space_dimension_, bool level_dof_access>
102class DoFCellAccessor;
103
104#endif
105
322template <int dim, int spacedim = dim>
325{
330
331public:
344 using cell_accessor = typename ActiveSelector::CellAccessor;
345
358 using face_accessor = typename ActiveSelector::FaceAccessor;
359
367 using line_iterator = typename ActiveSelector::line_iterator;
368
382 using active_line_iterator = typename ActiveSelector::active_line_iterator;
383
391 using quad_iterator = typename ActiveSelector::quad_iterator;
392
406 using active_quad_iterator = typename ActiveSelector::active_quad_iterator;
407
415 using hex_iterator = typename ActiveSelector::hex_iterator;
416
426 using active_hex_iterator = typename ActiveSelector::active_hex_iterator;
427
448 using active_cell_iterator = typename ActiveSelector::active_cell_iterator;
449
476 using cell_iterator = typename ActiveSelector::cell_iterator;
477
494 using face_iterator = typename ActiveSelector::face_iterator;
495
507 using active_face_iterator = typename ActiveSelector::active_face_iterator;
508
509 using level_cell_accessor = typename LevelSelector::CellAccessor;
510 using level_face_accessor = typename LevelSelector::FaceAccessor;
511
512 using level_cell_iterator = typename LevelSelector::cell_iterator;
513 using level_face_iterator = typename LevelSelector::face_iterator;
514
515
519 static constexpr unsigned int dimension = dim;
520
524 static constexpr unsigned int space_dimension = spacedim;
525
529 static const types::fe_index default_fe_index = 0;
530
534 using offset_type = unsigned int;
535
541
546
553 DoFHandler(const DoFHandler &) = delete;
554
558 virtual ~DoFHandler() override;
559
566 DoFHandler &
567 operator=(const DoFHandler &) = delete;
568
583 void
584 set_active_fe_indices(const std::vector<types::fe_index> &active_fe_indices);
585
592 void
593 set_active_fe_indices(const std::vector<unsigned int> &active_fe_indices);
594
608 std::vector<types::fe_index>
610
628 void
629 get_active_fe_indices(std::vector<unsigned int> &active_fe_indices) const;
630
644 void
645 set_future_fe_indices(const std::vector<types::fe_index> &future_fe_indices);
646
657 std::vector<types::fe_index>
659
667 void
669
699 void
701
705 void
707
713 void
715
719 bool
721
727 bool
729
744 bool
746
753 void
755
759 void
761
814 void
815 renumber_dofs(const std::vector<types::global_dof_index> &new_numbers);
816
821 void
822 renumber_dofs(const unsigned int level,
823 const std::vector<types::global_dof_index> &new_numbers);
824
853 unsigned int
855
868 unsigned int
870
871 /*--------------------------------------*/
872
883 begin(const unsigned int level = 0) const;
884
902 begin_active(const unsigned int level = 0) const;
903
909 end() const;
910
916 end(const unsigned int level) const;
917
924 end_active(const unsigned int level) const;
925
931 begin_mg(const unsigned int level = 0) const;
932
938 end_mg(const unsigned int level) const;
939
945 end_mg() const;
965
1008
1022
1039 cell_iterators_on_level(const unsigned int level) const;
1040
1057 active_cell_iterators_on_level(const unsigned int level) const;
1058
1075 mg_cell_iterators_on_level(const unsigned int level) const;
1076
1080 /*---------------------------------------*/
1081
1082
1107 n_dofs() const;
1108
1117 n_dofs(const unsigned int level) const;
1118
1125
1137 template <typename number>
1140 const std::map<types::boundary_id, const Function<spacedim, number> *>
1141 &boundary_ids) const;
1142
1148 n_boundary_dofs(const std::set<types::boundary_id> &boundary_ids) const;
1149
1165 const BlockInfo &
1166 block_info() const;
1167
1189
1195 const IndexSet &
1197
1202 const IndexSet &
1203 locally_owned_mg_dofs(const unsigned int level) const;
1204
1210 get_fe(const types::fe_index index = 0) const;
1211
1218
1224
1228 MPI_Comm
1230
1236 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
1237 "Access the MPI communicator with get_mpi_communicator() instead.")
1238 MPI_Comm
1239 get_communicator() const;
1240
1257 void
1258 prepare_for_serialization_of_active_fe_indices();
1259
1275 void
1276 deserialize_active_fe_indices();
1277
1286 virtual std::size_t
1287 memory_consumption() const;
1288
1294 template <class Archive>
1295 void
1296 save(Archive &ar, const unsigned int version) const;
1297
1303 template <class Archive>
1304 void
1305 load(Archive &ar, const unsigned int version);
1306
1307#ifdef DOXYGEN
1313 template <class Archive>
1314 void
1315 serialize(Archive &archive, const unsigned int version);
1316#else
1317 // This macro defines the serialize() method that is compatible with
1318 // the templated save() and load() method that have been implemented.
1319 BOOST_SERIALIZATION_SPLIT_MEMBER()
1320#endif
1321
1325 DeclException0(ExcNoFESelected);
1330 DeclException0(ExcInvalidBoundaryIndicator);
1335 DeclException1(ExcInvalidLevel,
1336 int,
1337 << "The given level " << arg1
1338 << " is not in the valid range!");
1343 DeclException1(ExcNewNumbersNotConsecutive,
1345 << "The given list of new dof indices is not consecutive: "
1346 << "the index " << arg1 << " does not exist.");
1350 DeclException2(ExcInvalidFEIndex,
1351 int,
1352 int,
1353 << "The mesh contains a cell with an active FE index of "
1354 << arg1 << ", but the finite element collection only has "
1355 << arg2 << " elements");
1356
1361 DeclExceptionMsg(ExcOnlyAvailableWithHP,
1362 "The current function doesn't make sense when used with a "
1363 "DoFHandler without hp-capabilities.");
1364
1369 DeclExceptionMsg(ExcNotImplementedWithHP,
1370 "The current function has not yet been implemented for a "
1371 "DoFHandler with hp-capabilities.");
1372
1373private:
1381 {
1382 public:
1387
1393 void
1394 init(const unsigned int coarsest_level,
1395 const unsigned int finest_level,
1396 const unsigned int dofs_per_vertex);
1397
1401 unsigned int
1403
1407 unsigned int
1409
1415 access_index(const unsigned int level,
1416 const unsigned int dof_number,
1417 const unsigned int dofs_per_vertex);
1418
1419 private:
1423 unsigned int coarsest_level;
1424
1428 unsigned int finest_level;
1429
1439 std::unique_ptr<types::global_dof_index[]> indices;
1440 };
1441
1449 {
1454 std::map<const cell_iterator, const types::fe_index>
1456
1461 std::map<const cell_iterator, const types::fe_index> refined_cells_fe_index;
1462
1467 std::map<const cell_iterator, const types::fe_index>
1469
1475 std::vector<types::fe_index> active_fe_indices;
1476
1482 std::unique_ptr<
1483 parallel::distributed::
1484 CellDataTransfer<dim, spacedim, std::vector<types::fe_index>>>
1486 };
1487
1492
1498
1504
1511
1516 std::unique_ptr<::internal::DoFHandlerImplementation::Policy::
1517 PolicyBase<dim, spacedim>>
1519
1528
1532 std::vector<::internal::DoFHandlerImplementation::NumberCache>
1534
1540 mutable std::vector<std::array<std::vector<types::global_dof_index>, dim + 1>>
1542
1551 mutable std::vector<std::array<std::vector<offset_type>, dim + 1>>
1553
1559 mutable std::array<std::vector<types::fe_index>, dim + 1>
1561
1565 mutable std::array<std::vector<offset_type>, dim + 1> hp_object_fe_ptr;
1566
1571 mutable std::vector<std::vector<types::fe_index>> hp_cell_active_fe_indices;
1572
1577 mutable std::vector<std::vector<types::fe_index>> hp_cell_future_fe_indices;
1578
1583 std::vector<MGVertexDoFs> mg_vertex_dofs;
1584
1588 std::vector<
1589 std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel<dim>>>
1591
1595 std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces<dim>>
1597
1602 std::unique_ptr<ActiveFEIndexTransfer> active_fe_index_transfer;
1603
1608 std::vector<boost::signals2::connection> tria_listeners;
1609
1615 std::vector<boost::signals2::connection> tria_listeners_for_transfer;
1616
1620 void
1622
1626 void
1628
1632 void
1634
1638 void
1640
1652 void
1654
1667 void
1669
1679 void
1681
1690 void
1692
1701 void
1703
1712 void
1714
1715
1716 // Make accessor objects friends.
1717 template <int, int, int, bool>
1718 friend class ::DoFAccessor;
1719 template <int, int, bool>
1720 friend class ::DoFCellAccessor;
1721 friend struct ::internal::DoFAccessorImplementation::Implementation;
1722 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1723
1724 // Likewise for DoFLevel objects since they need to access the vertex dofs
1725 // in the functions that set and retrieve vertex dof indices.
1726 friend struct ::internal::DoFHandlerImplementation::Implementation;
1727 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1728 friend struct ::internal::DoFHandlerImplementation::Policy::
1729 Implementation;
1730
1731 // explicitly check for sensible template arguments, but not on windows
1732 // because MSVC creates bogus warnings during normal compilation
1733#ifndef DEAL_II_MSVC
1734 static_assert(dim <= spacedim,
1735 "The dimension <dim> of a DoFHandler must be less than or "
1736 "equal to the space dimension <spacedim> in which it lives.");
1737#endif
1738};
1739
1740namespace internal
1741{
1742 namespace hp
1743 {
1744 namespace DoFHandlerImplementation
1745 {
1757 template <int dim, int spacedim>
1758 void
1760
1785 template <int dim, int spacedim = dim>
1786 unsigned int
1788 const typename DoFHandler<dim, spacedim>::cell_iterator &parent);
1789
1795 "No FiniteElement has been found in your FECollection that is "
1796 "dominated by all children of a cell you are trying to coarsen!");
1797 } // namespace DoFHandlerImplementation
1798 } // namespace hp
1799} // namespace internal
1800
1801#ifndef DOXYGEN
1802
1803/* ----------------------- Inline functions ----------------------------------
1804 */
1805
1806
1807template <int dim, int spacedim>
1810{
1811 return hp_capability_enabled;
1812}
1813
1814
1815
1816template <int dim, int spacedim>
1819{
1820 return mg_number_cache.size() > 0;
1821}
1822
1823
1824
1825template <int dim, int spacedim>
1828{
1829 return number_cache.n_global_dofs > 0;
1830}
1831
1832
1833
1834template <int dim, int spacedim>
1837{
1838 return number_cache.n_global_dofs;
1839}
1840
1841
1842
1843template <int dim, int spacedim>
1846 DoFHandler<dim, spacedim>::n_dofs(const unsigned int level) const
1847{
1848 Assert(has_level_dofs(),
1849 ExcMessage(
1850 "n_dofs(level) can only be called after distribute_mg_dofs()"));
1851 Assert(level < mg_number_cache.size(), ExcInvalidLevel(level));
1852 return mg_number_cache[level].n_global_dofs;
1853}
1854
1855
1856
1857template <int dim, int spacedim>
1860{
1861 return number_cache.n_locally_owned_dofs;
1862}
1863
1864
1865
1866template <int dim, int spacedim>
1869{
1870 return number_cache.locally_owned_dofs;
1871}
1872
1873
1874
1875template <int dim, int spacedim>
1878 const unsigned int level) const
1879{
1880 Assert(level < this->get_triangulation().n_global_levels(),
1881 ExcMessage("The given level index exceeds the number of levels "
1882 "present in the triangulation"));
1883 Assert(
1884 mg_number_cache.size() == this->get_triangulation().n_global_levels(),
1885 ExcMessage(
1886 "The level dofs are not set up properly! Did you call distribute_mg_dofs()?"));
1887 return mg_number_cache[level].locally_owned_dofs;
1888}
1889
1890
1891
1892template <int dim, int spacedim>
1895 const types::fe_index number) const
1896{
1897 Assert(fe_collection.size() > 0,
1898 ExcMessage("No finite element collection is associated with "
1899 "this DoFHandler"));
1900 return fe_collection[number];
1901}
1902
1903
1904
1905template <int dim, int spacedim>
1909{
1910 return fe_collection;
1911}
1912
1913
1914
1915template <int dim, int spacedim>
1919{
1920 Assert(tria != nullptr,
1921 ExcMessage("This DoFHandler object has not been associated "
1922 "with a triangulation."));
1923 return *tria;
1924}
1925
1926
1927
1928template <int dim, int spacedim>
1931{
1932 Assert(tria != nullptr,
1933 ExcMessage("This DoFHandler object has not been associated "
1934 "with a triangulation."));
1935 return tria->get_mpi_communicator();
1936}
1937
1938
1939
1940template <int dim, int spacedim>
1943{
1944 return get_mpi_communicator();
1945}
1946
1947
1948
1949template <int dim, int spacedim>
1952{
1953 Assert(this->hp_capability_enabled == false, ExcNotImplementedWithHP());
1954
1955 return block_info_object;
1956}
1957
1958
1959
1960template <int dim, int spacedim>
1962template <typename number>
1964 const std::map<types::boundary_id, const Function<spacedim, number> *>
1965 &boundary_ids) const
1966{
1967 Assert(!(dim == 2 && spacedim == 3) || this->hp_capability_enabled == false,
1968 ExcNotImplementedWithHP());
1969
1970 // extract the set of boundary ids and forget about the function object
1971 // pointers
1972 std::set<types::boundary_id> boundary_ids_only;
1973 for (typename std::map<types::boundary_id,
1974 const Function<spacedim, number> *>::const_iterator p =
1975 boundary_ids.begin();
1976 p != boundary_ids.end();
1977 ++p)
1978 boundary_ids_only.insert(p->first);
1979
1980 // then just hand everything over to the other function that does the work
1981 return n_boundary_dofs(boundary_ids_only);
1982}
1983
1984
1985
1986namespace internal
1987{
1995 template <int dim, int spacedim>
1996 std::string
1997 policy_to_string(const ::internal::DoFHandlerImplementation::Policy::
1998 PolicyBase<dim, spacedim> &policy);
1999} // namespace internal
2000
2001
2002
2003template <int dim, int spacedim>
2005template <class Archive>
2006void DoFHandler<dim, spacedim>::save(Archive &ar, const unsigned int) const
2007{
2008 if (this->hp_capability_enabled)
2009 {
2010 ar &this->object_dof_indices;
2011 ar &this->object_dof_ptr;
2012
2013 ar &this->hp_cell_active_fe_indices;
2014 ar &this->hp_cell_future_fe_indices;
2015
2016 ar &hp_object_fe_ptr;
2017 ar &hp_object_fe_indices;
2018
2019 ar &number_cache;
2020
2021 ar &mg_number_cache;
2022
2023 // write out the number of triangulation cells and later check during
2024 // loading that this number is indeed correct; same with something that
2025 // identifies the policy
2026 const unsigned int n_cells = this->tria->n_cells();
2027 std::string policy_name =
2028 ::internal::policy_to_string(*this->policy);
2029
2030 ar &n_cells &policy_name;
2031 }
2032 else
2033 {
2034 ar &this->block_info_object;
2035 ar &number_cache;
2036
2037 ar &this->object_dof_indices;
2038 ar &this->object_dof_ptr;
2039
2040 // write out the number of triangulation cells and later check during
2041 // loading that this number is indeed correct; same with something that
2042 // identifies the FE and the policy
2043 unsigned int n_cells = this->tria->n_cells();
2044 std::string fe_name = this->get_fe(0).get_name();
2045 std::string policy_name = internal::policy_to_string(*this->policy);
2046
2047 ar &n_cells &fe_name &policy_name;
2048 }
2049}
2050
2051
2052
2053template <int dim, int spacedim>
2055template <class Archive>
2056void DoFHandler<dim, spacedim>::load(Archive &ar, const unsigned int)
2057{
2058 if (this->hp_capability_enabled)
2059 {
2060 ar &this->object_dof_indices;
2061 ar &this->object_dof_ptr;
2062
2063 ar &this->hp_cell_active_fe_indices;
2064 ar &this->hp_cell_future_fe_indices;
2065
2066 ar &hp_object_fe_ptr;
2067 ar &hp_object_fe_indices;
2068
2069 ar &number_cache;
2070
2071 ar &mg_number_cache;
2072
2073 // these are the checks that correspond to the last block in the save()
2074 // function
2075 unsigned int n_cells;
2076 std::string policy_name;
2077
2078 ar &n_cells &policy_name;
2079
2081 n_cells == this->tria->n_cells(),
2082 ExcMessage(
2083 "The object being loaded into does not match the triangulation "
2084 "that has been stored previously."));
2086 policy_name == ::internal::policy_to_string(*this->policy),
2087 ExcMessage("The policy currently associated with this DoFHandler (" +
2088 ::internal::policy_to_string(*this->policy) +
2089 ") does not match the one that was associated with the "
2090 "DoFHandler previously stored (" +
2091 policy_name + ")."));
2092 }
2093 else
2094 {
2095 ar &this->block_info_object;
2096 ar &number_cache;
2097
2098 object_dof_indices.clear();
2099
2100 object_dof_ptr.clear();
2101
2102 ar &this->object_dof_indices;
2103 ar &this->object_dof_ptr;
2104
2105 // these are the checks that correspond to the last block in the save()
2106 // function
2107 unsigned int n_cells;
2108 std::string fe_name;
2109 std::string policy_name;
2110
2111 ar &n_cells &fe_name &policy_name;
2112
2114 n_cells == this->tria->n_cells(),
2115 ExcMessage(
2116 "The object being loaded into does not match the triangulation "
2117 "that has been stored previously."));
2119 fe_name == this->get_fe(0).get_name(),
2120 ExcMessage(
2121 "The finite element associated with this DoFHandler does not match "
2122 "the one that was associated with the DoFHandler previously stored."));
2123 AssertThrow(policy_name == internal::policy_to_string(*this->policy),
2124 ExcMessage(
2125 "The policy currently associated with this DoFHandler (" +
2126 internal::policy_to_string(*this->policy) +
2127 ") does not match the one that was associated with the "
2128 "DoFHandler previously stored (" +
2129 policy_name + ")."));
2130 }
2131}
2132
2133
2134
2135template <int dim, int spacedim>
2139 const unsigned int level,
2140 const unsigned int dof_number,
2141 const unsigned int dofs_per_vertex)
2142{
2145 return indices[dofs_per_vertex * (level - coarsest_level) + dof_number];
2146}
2147
2148
2149
2150extern template class DoFHandler<1, 1>;
2151extern template class DoFHandler<1, 2>;
2152extern template class DoFHandler<1, 3>;
2153extern template class DoFHandler<2, 2>;
2154extern template class DoFHandler<2, 3>;
2155extern template class DoFHandler<3, 3>;
2156
2157
2158#endif // DOXYGEN
2159
2161
2162#endif
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition block_info.h:97
static const unsigned int dim
static const unsigned int spacedim
types::global_dof_index & access_index(const unsigned int level, const unsigned int dof_number, const unsigned int dofs_per_vertex)
void init(const unsigned int coarsest_level, const unsigned int finest_level, const unsigned int dofs_per_vertex)
unsigned int get_finest_level() const
std::unique_ptr< types::global_dof_index[]> indices
unsigned int get_coarsest_level() const
cell_iterator end() const
void pre_distributed_transfer_action()
void post_transfer_action()
std::vector< types::fe_index > get_future_fe_indices() const
unsigned int max_couplings_between_dofs() const
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > mg_levels
std::unique_ptr< ActiveFEIndexTransfer > active_fe_index_transfer
hp::FECollection< dim, spacedim > fe_collection
void create_active_fe_table()
std::vector< std::vector< types::fe_index > > hp_cell_future_fe_indices
const hp::FECollection< dim, spacedim > & get_fe_collection() const
void save(Archive &ar, const unsigned int version) const
void set_active_fe_indices(const std::vector< unsigned int > &active_fe_indices)
BlockInfo block_info_object
bool has_active_dofs() const
std::vector< boost::signals2::connection > tria_listeners_for_transfer
level_cell_iterator end_mg() const
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
void renumber_dofs(const unsigned int level, const std::vector< types::global_dof_index > &new_numbers)
level_cell_iterator begin_mg(const unsigned int level=0) const
void distribute_dofs(const FiniteElement< dim, spacedim > &fe)
types::global_dof_index n_boundary_dofs() const
DoFHandler(const Triangulation< dim, spacedim > &tria)
void clear_mg_space()
void connect_to_triangulation_signals()
std::vector< std::vector< types::fe_index > > hp_cell_active_fe_indices
std::vector< MGVertexDoFs > mg_vertex_dofs
std::vector< boost::signals2::connection > tria_listeners
void post_distributed_transfer_action()
std::vector< std::array< std::vector< types::global_dof_index >, dim+1 > > object_dof_indices
const Triangulation< dim, spacedim > & get_triangulation() const
typename LevelSelector::CellAccessor level_cell_accessor
void pre_transfer_action()
void clear_space()
const IndexSet & locally_owned_dofs() const
void reinit(const Triangulation< dim, spacedim > &tria)
std::vector< std::array< std::vector< offset_type >, dim+1 > > object_dof_ptr
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > mg_faces
void serialize(Archive &archive, const unsigned int version)
DoFHandler & operator=(const DoFHandler &)=delete
DoFHandler(const DoFHandler &)=delete
active_cell_iterator begin_active(const unsigned int level=0) const
std::vector< types::fe_index > get_active_fe_indices() const
void distribute_mg_dofs()
active_cell_iterator end_active(const unsigned int level) const
bool hp_capability_enabled
bool has_hp_capabilities() const
void initialize_local_block_info()
types::global_dof_index n_dofs() const
void update_active_fe_table()
level_cell_iterator end_mg(const unsigned int level) const
cell_iterator end(const unsigned int level) const
typename LevelSelector::cell_iterator level_cell_iterator
std::vector<::internal::DoFHandlerImplementation::NumberCache > mg_number_cache
void load(Archive &ar, const unsigned int version)
typename LevelSelector::face_iterator level_face_iterator
cell_iterator begin(const unsigned int level=0) const
MPI_Comm get_mpi_communicator() const
std::array< std::vector< offset_type >, dim+1 > hp_object_fe_ptr
MPI_Comm get_communicator() const
void clear()
void get_active_fe_indices(std::vector< unsigned int > &active_fe_indices) const
void distribute_dofs(const hp::FECollection< dim, spacedim > &fe)
types::global_dof_index n_boundary_dofs(const std::map< types::boundary_id, const Function< spacedim, number > * > &boundary_ids) const
ObserverPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
bool has_level_dofs() const
std::unique_ptr<::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > > policy
typename LevelSelector::FaceAccessor level_face_accessor
types::global_dof_index n_dofs(const unsigned int level) const
std::array< std::vector< types::fe_index >, dim+1 > hp_object_fe_indices
void set_active_fe_indices(const std::vector< types::fe_index > &active_fe_indices)
unsigned int max_couplings_between_boundary_dofs() const
void setup_policy()
types::global_dof_index n_locally_owned_dofs() const
types::global_dof_index n_boundary_dofs(const std::set< types::boundary_id > &boundary_ids) const
const BlockInfo & block_info() const
::internal::DoFHandlerImplementation::NumberCache number_cache
void set_future_fe_indices(const std::vector< types::fe_index > &future_fe_indices)
virtual ~DoFHandler() override
#define DEAL_II_DEPRECATED
Definition config.h:280
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:242
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
unsigned int level
Definition grid_out.cc:4635
IteratorRange< active_cell_iterator > active_cell_iterators() const
IteratorRange< level_cell_iterator > mg_cell_iterators() const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
IteratorRange< cell_iterator > cell_iterators() const
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
#define DeclException0(Exception0)
#define Assert(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcNoDominatedFiniteElementOnChildren()
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcInvalidLevel(int arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
typename ActiveSelector::FaceAccessor face_accessor
typename ActiveSelector::cell_iterator cell_iterator
typename ActiveSelector::line_iterator line_iterator
typename ActiveSelector::face_iterator face_iterator
typename ActiveSelector::active_line_iterator active_line_iterator
typename ActiveSelector::active_cell_iterator active_cell_iterator
typename ActiveSelector::CellAccessor cell_accessor
typename ActiveSelector::quad_iterator quad_iterator
typename ActiveSelector::active_hex_iterator active_hex_iterator
typename ActiveSelector::active_quad_iterator active_quad_iterator
typename ActiveSelector::hex_iterator hex_iterator
typename ActiveSelector::active_face_iterator active_face_iterator
Definition hp.h:117
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition tria.cc:14904
unsigned int dominated_future_fe_on_children(const typename DoFHandler< dim, spacedim >::cell_iterator &parent)
void communicate_future_fe_indices(DoFHandler< dim, spacedim > &dof_handler)
std::string policy_to_string(const ::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > &policy)
STL namespace.
unsigned short int fe_index
Definition types.h:72
std::map< const cell_iterator, const types::fe_index > refined_cells_fe_index
std::map< const cell_iterator, const types::fe_index > persisting_cells_fe_index
std::vector< types::fe_index > active_fe_indices
std::map< const cell_iterator, const types::fe_index > coarsened_cells_fe_index
std::unique_ptr< parallel::distributed::CellDataTransfer< dim, spacedim, std::vector< types::fe_index > > > cell_data_transfer