deal.II version GIT relicensing-2017-g7ae82fad8b 2024-10-22 19:20: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_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
40#include <map>
41#include <memory>
42#include <set>
43#include <vector>
44
46
47// Forward declarations
48#ifndef DOXYGEN
49template <int dim, int spacedim>
50class FiniteElement;
51template <int dim, int spacedim>
53class Triangulation;
54
55namespace internal
56{
57 namespace DoFHandlerImplementation
58 {
59 struct Implementation;
60
61 namespace Policy
62 {
63 template <int dim, int spacedim>
64 class PolicyBase;
65 struct Implementation;
66 } // namespace Policy
67 } // namespace DoFHandlerImplementation
68
69 namespace DoFAccessorImplementation
70 {
71 struct Implementation;
72 }
73
74 namespace DoFCellAccessorImplementation
75 {
76 struct Implementation;
77 }
78
79 namespace hp
80 {
81 namespace DoFHandlerImplementation
82 {
83 struct Implementation;
84 }
85 } // namespace hp
86} // namespace internal
87
88namespace parallel
89{
90 namespace distributed
91 {
92 template <int dim, int spacedim, typename VectorType>
93 class CellDataTransfer;
94 }
95} // namespace parallel
96#endif
97
314template <int dim, int spacedim = dim>
317{
322
323public:
336 using cell_accessor = typename ActiveSelector::CellAccessor;
337
350 using face_accessor = typename ActiveSelector::FaceAccessor;
351
359 using line_iterator = typename ActiveSelector::line_iterator;
360
374 using active_line_iterator = typename ActiveSelector::active_line_iterator;
375
383 using quad_iterator = typename ActiveSelector::quad_iterator;
384
398 using active_quad_iterator = typename ActiveSelector::active_quad_iterator;
399
407 using hex_iterator = typename ActiveSelector::hex_iterator;
408
418 using active_hex_iterator = typename ActiveSelector::active_hex_iterator;
419
440 using active_cell_iterator = typename ActiveSelector::active_cell_iterator;
441
468 using cell_iterator = typename ActiveSelector::cell_iterator;
469
486 using face_iterator = typename ActiveSelector::face_iterator;
487
499 using active_face_iterator = typename ActiveSelector::active_face_iterator;
500
501 using level_cell_accessor = typename LevelSelector::CellAccessor;
502 using level_face_accessor = typename LevelSelector::FaceAccessor;
503
504 using level_cell_iterator = typename LevelSelector::cell_iterator;
505 using level_face_iterator = typename LevelSelector::face_iterator;
506
507
511 static constexpr unsigned int dimension = dim;
512
516 static constexpr unsigned int space_dimension = spacedim;
517
521 static const types::fe_index default_fe_index = 0;
522
526 using offset_type = unsigned int;
527
533
538
545 DoFHandler(const DoFHandler &) = delete;
546
550 virtual ~DoFHandler() override;
551
558 DoFHandler &
559 operator=(const DoFHandler &) = delete;
560
575 void
576 set_active_fe_indices(const std::vector<types::fe_index> &active_fe_indices);
577
584 void
585 set_active_fe_indices(const std::vector<unsigned int> &active_fe_indices);
586
600 std::vector<types::fe_index>
602
620 void
621 get_active_fe_indices(std::vector<unsigned int> &active_fe_indices) const;
622
636 void
637 set_future_fe_indices(const std::vector<types::fe_index> &future_fe_indices);
638
649 std::vector<types::fe_index>
651
659 void
661
691 void
693
697 void
699
705 void
707
711 bool
713
719 bool
721
736 bool
738
745 void
747
751 void
753
806 void
807 renumber_dofs(const std::vector<types::global_dof_index> &new_numbers);
808
813 void
814 renumber_dofs(const unsigned int level,
815 const std::vector<types::global_dof_index> &new_numbers);
816
845 unsigned int
847
860 unsigned int
862
863 /*--------------------------------------*/
864
875 begin(const unsigned int level = 0) const;
876
894 begin_active(const unsigned int level = 0) const;
895
901 end() const;
902
908 end(const unsigned int level) const;
909
916 end_active(const unsigned int level) const;
917
923 begin_mg(const unsigned int level = 0) const;
924
930 end_mg(const unsigned int level) const;
931
937 end_mg() const;
957
1000
1014
1031 cell_iterators_on_level(const unsigned int level) const;
1032
1049 active_cell_iterators_on_level(const unsigned int level) const;
1050
1067 mg_cell_iterators_on_level(const unsigned int level) const;
1068
1072 /*---------------------------------------*/
1073
1074
1099 n_dofs() const;
1100
1109 n_dofs(const unsigned int level) const;
1110
1117
1129 template <typename number>
1132 const std::map<types::boundary_id, const Function<spacedim, number> *>
1133 &boundary_ids) const;
1134
1140 n_boundary_dofs(const std::set<types::boundary_id> &boundary_ids) const;
1141
1157 const BlockInfo &
1158 block_info() const;
1159
1181
1187 const IndexSet &
1189
1194 const IndexSet &
1195 locally_owned_mg_dofs(const unsigned int level) const;
1196
1202 get_fe(const types::fe_index index = 0) const;
1203
1210
1216
1220 MPI_Comm
1222
1228 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
1229 "Access the MPI communicator with get_mpi_communicator() instead.")
1230 MPI_Comm
1231 get_communicator() const;
1232
1249 void
1250 prepare_for_serialization_of_active_fe_indices();
1251
1267 void
1268 deserialize_active_fe_indices();
1269
1278 virtual std::size_t
1279 memory_consumption() const;
1280
1286 template <class Archive>
1287 void
1288 save(Archive &ar, const unsigned int version) const;
1289
1295 template <class Archive>
1296 void
1297 load(Archive &ar, const unsigned int version);
1298
1299#ifdef DOXYGEN
1305 template <class Archive>
1306 void
1307 serialize(Archive &archive, const unsigned int version);
1308#else
1309 // This macro defines the serialize() method that is compatible with
1310 // the templated save() and load() method that have been implemented.
1311 BOOST_SERIALIZATION_SPLIT_MEMBER()
1312#endif
1313
1317 DeclException0(ExcNoFESelected);
1322 DeclException0(ExcInvalidBoundaryIndicator);
1327 DeclException1(ExcInvalidLevel,
1328 int,
1329 << "The given level " << arg1
1330 << " is not in the valid range!");
1335 DeclException1(ExcNewNumbersNotConsecutive,
1337 << "The given list of new dof indices is not consecutive: "
1338 << "the index " << arg1 << " does not exist.");
1342 DeclException2(ExcInvalidFEIndex,
1343 int,
1344 int,
1345 << "The mesh contains a cell with an active FE index of "
1346 << arg1 << ", but the finite element collection only has "
1347 << arg2 << " elements");
1348
1353 DeclExceptionMsg(ExcOnlyAvailableWithHP,
1354 "The current function doesn't make sense when used with a "
1355 "DoFHandler without hp-capabilities.");
1356
1361 DeclExceptionMsg(ExcNotImplementedWithHP,
1362 "The current function has not yet been implemented for a "
1363 "DoFHandler with hp-capabilities.");
1364
1365private:
1373 {
1374 public:
1379
1385 void
1386 init(const unsigned int coarsest_level,
1387 const unsigned int finest_level,
1388 const unsigned int dofs_per_vertex);
1389
1393 unsigned int
1395
1399 unsigned int
1401
1407 access_index(const unsigned int level,
1408 const unsigned int dof_number,
1409 const unsigned int dofs_per_vertex);
1410
1411 private:
1415 unsigned int coarsest_level;
1416
1420 unsigned int finest_level;
1421
1431 std::unique_ptr<types::global_dof_index[]> indices;
1432 };
1433
1441 {
1446 std::map<const cell_iterator, const types::fe_index>
1448
1453 std::map<const cell_iterator, const types::fe_index> refined_cells_fe_index;
1454
1459 std::map<const cell_iterator, const types::fe_index>
1461
1467 std::vector<types::fe_index> active_fe_indices;
1468
1474 std::unique_ptr<
1475 parallel::distributed::
1476 CellDataTransfer<dim, spacedim, std::vector<types::fe_index>>>
1478 };
1479
1484
1490
1496
1503
1508 std::unique_ptr<::internal::DoFHandlerImplementation::Policy::
1509 PolicyBase<dim, spacedim>>
1511
1520
1524 std::vector<::internal::DoFHandlerImplementation::NumberCache>
1526
1532 mutable std::vector<std::array<std::vector<types::global_dof_index>, dim + 1>>
1534
1543 mutable std::vector<std::array<std::vector<offset_type>, dim + 1>>
1545
1551 mutable std::array<std::vector<types::fe_index>, dim + 1>
1553
1557 mutable std::array<std::vector<offset_type>, dim + 1> hp_object_fe_ptr;
1558
1563 mutable std::vector<std::vector<types::fe_index>> hp_cell_active_fe_indices;
1564
1569 mutable std::vector<std::vector<types::fe_index>> hp_cell_future_fe_indices;
1570
1575 std::vector<MGVertexDoFs> mg_vertex_dofs;
1576
1580 std::vector<
1581 std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel<dim>>>
1583
1587 std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces<dim>>
1589
1594 std::unique_ptr<ActiveFEIndexTransfer> active_fe_index_transfer;
1595
1600 std::vector<boost::signals2::connection> tria_listeners;
1601
1607 std::vector<boost::signals2::connection> tria_listeners_for_transfer;
1608
1612 void
1614
1618 void
1620
1624 void
1626
1630 void
1632
1644 void
1646
1659 void
1661
1671 void
1673
1682 void
1684
1693 void
1695
1704 void
1706
1707
1708 // Make accessor objects friends.
1709 template <int, int, int, bool>
1710 friend class ::DoFAccessor;
1711 template <int, int, bool>
1712 friend class ::DoFCellAccessor;
1713 friend struct ::internal::DoFAccessorImplementation::Implementation;
1714 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1715
1716 // Likewise for DoFLevel objects since they need to access the vertex dofs
1717 // in the functions that set and retrieve vertex dof indices.
1718 friend struct ::internal::DoFHandlerImplementation::Implementation;
1719 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1720 friend struct ::internal::DoFHandlerImplementation::Policy::
1721 Implementation;
1722
1723 // explicitly check for sensible template arguments, but not on windows
1724 // because MSVC creates bogus warnings during normal compilation
1725#ifndef DEAL_II_MSVC
1726 static_assert(dim <= spacedim,
1727 "The dimension <dim> of a DoFHandler must be less than or "
1728 "equal to the space dimension <spacedim> in which it lives.");
1729#endif
1730};
1731
1732namespace internal
1733{
1734 namespace hp
1735 {
1736 namespace DoFHandlerImplementation
1737 {
1749 template <int dim, int spacedim>
1750 void
1752
1777 template <int dim, int spacedim = dim>
1778 unsigned int
1780 const typename DoFHandler<dim, spacedim>::cell_iterator &parent);
1781
1787 "No FiniteElement has been found in your FECollection that is "
1788 "dominated by all children of a cell you are trying to coarsen!");
1789 } // namespace DoFHandlerImplementation
1790 } // namespace hp
1791} // namespace internal
1792
1793#ifndef DOXYGEN
1794
1795/* ----------------------- Inline functions ----------------------------------
1796 */
1797
1798
1799template <int dim, int spacedim>
1802{
1803 return hp_capability_enabled;
1804}
1805
1806
1807
1808template <int dim, int spacedim>
1811{
1812 return mg_number_cache.size() > 0;
1813}
1814
1815
1816
1817template <int dim, int spacedim>
1820{
1821 return number_cache.n_global_dofs > 0;
1822}
1823
1824
1825
1826template <int dim, int spacedim>
1829{
1830 return number_cache.n_global_dofs;
1831}
1832
1833
1834
1835template <int dim, int spacedim>
1838 DoFHandler<dim, spacedim>::n_dofs(const unsigned int level) const
1839{
1840 Assert(has_level_dofs(),
1841 ExcMessage(
1842 "n_dofs(level) can only be called after distribute_mg_dofs()"));
1843 Assert(level < mg_number_cache.size(), ExcInvalidLevel(level));
1844 return mg_number_cache[level].n_global_dofs;
1845}
1846
1847
1848
1849template <int dim, int spacedim>
1852{
1853 return number_cache.n_locally_owned_dofs;
1854}
1855
1856
1857
1858template <int dim, int spacedim>
1861{
1862 return number_cache.locally_owned_dofs;
1863}
1864
1865
1866
1867template <int dim, int spacedim>
1870 const unsigned int level) const
1871{
1872 Assert(level < this->get_triangulation().n_global_levels(),
1873 ExcMessage("The given level index exceeds the number of levels "
1874 "present in the triangulation"));
1875 Assert(
1876 mg_number_cache.size() == this->get_triangulation().n_global_levels(),
1877 ExcMessage(
1878 "The level dofs are not set up properly! Did you call distribute_mg_dofs()?"));
1879 return mg_number_cache[level].locally_owned_dofs;
1880}
1881
1882
1883
1884template <int dim, int spacedim>
1887 const types::fe_index number) const
1888{
1889 Assert(fe_collection.size() > 0,
1890 ExcMessage("No finite element collection is associated with "
1891 "this DoFHandler"));
1892 return fe_collection[number];
1893}
1894
1895
1896
1897template <int dim, int spacedim>
1901{
1902 return fe_collection;
1903}
1904
1905
1906
1907template <int dim, int spacedim>
1911{
1912 Assert(tria != nullptr,
1913 ExcMessage("This DoFHandler object has not been associated "
1914 "with a triangulation."));
1915 return *tria;
1916}
1917
1918
1919
1920template <int dim, int spacedim>
1923{
1924 Assert(tria != nullptr,
1925 ExcMessage("This DoFHandler object has not been associated "
1926 "with a triangulation."));
1927 return tria->get_mpi_communicator();
1928}
1929
1930
1931
1932template <int dim, int spacedim>
1935{
1936 return get_mpi_communicator();
1937}
1938
1939
1940
1941template <int dim, int spacedim>
1944{
1945 Assert(this->hp_capability_enabled == false, ExcNotImplementedWithHP());
1946
1947 return block_info_object;
1948}
1949
1950
1951
1952template <int dim, int spacedim>
1954template <typename number>
1956 const std::map<types::boundary_id, const Function<spacedim, number> *>
1957 &boundary_ids) const
1958{
1959 Assert(!(dim == 2 && spacedim == 3) || this->hp_capability_enabled == false,
1960 ExcNotImplementedWithHP());
1961
1962 // extract the set of boundary ids and forget about the function object
1963 // pointers
1964 std::set<types::boundary_id> boundary_ids_only;
1965 for (typename std::map<types::boundary_id,
1966 const Function<spacedim, number> *>::const_iterator p =
1967 boundary_ids.begin();
1968 p != boundary_ids.end();
1969 ++p)
1970 boundary_ids_only.insert(p->first);
1971
1972 // then just hand everything over to the other function that does the work
1973 return n_boundary_dofs(boundary_ids_only);
1974}
1975
1976
1977
1978namespace internal
1979{
1987 template <int dim, int spacedim>
1988 std::string
1989 policy_to_string(const ::internal::DoFHandlerImplementation::Policy::
1990 PolicyBase<dim, spacedim> &policy);
1991} // namespace internal
1992
1993
1994
1995template <int dim, int spacedim>
1997template <class Archive>
1998void DoFHandler<dim, spacedim>::save(Archive &ar, const unsigned int) const
1999{
2000 if (this->hp_capability_enabled)
2001 {
2002 ar &this->object_dof_indices;
2003 ar &this->object_dof_ptr;
2004
2005 ar &this->hp_cell_active_fe_indices;
2006 ar &this->hp_cell_future_fe_indices;
2007
2008 ar &hp_object_fe_ptr;
2009 ar &hp_object_fe_indices;
2010
2011 ar &number_cache;
2012
2013 ar &mg_number_cache;
2014
2015 // write out the number of triangulation cells and later check during
2016 // loading that this number is indeed correct; same with something that
2017 // identifies the policy
2018 const unsigned int n_cells = this->tria->n_cells();
2019 std::string policy_name =
2020 ::internal::policy_to_string(*this->policy);
2021
2022 ar &n_cells &policy_name;
2023 }
2024 else
2025 {
2026 ar &this->block_info_object;
2027 ar &number_cache;
2028
2029 ar &this->object_dof_indices;
2030 ar &this->object_dof_ptr;
2031
2032 // write out the number of triangulation cells and later check during
2033 // loading that this number is indeed correct; same with something that
2034 // identifies the FE and the policy
2035 unsigned int n_cells = this->tria->n_cells();
2036 std::string fe_name = this->get_fe(0).get_name();
2037 std::string policy_name = internal::policy_to_string(*this->policy);
2038
2039 ar &n_cells &fe_name &policy_name;
2040 }
2041}
2042
2043
2044
2045template <int dim, int spacedim>
2047template <class Archive>
2048void DoFHandler<dim, spacedim>::load(Archive &ar, const unsigned int)
2049{
2050 if (this->hp_capability_enabled)
2051 {
2052 ar &this->object_dof_indices;
2053 ar &this->object_dof_ptr;
2054
2055 ar &this->hp_cell_active_fe_indices;
2056 ar &this->hp_cell_future_fe_indices;
2057
2058 ar &hp_object_fe_ptr;
2059 ar &hp_object_fe_indices;
2060
2061 ar &number_cache;
2062
2063 ar &mg_number_cache;
2064
2065 // these are the checks that correspond to the last block in the save()
2066 // function
2067 unsigned int n_cells;
2068 std::string policy_name;
2069
2070 ar &n_cells &policy_name;
2071
2073 n_cells == this->tria->n_cells(),
2074 ExcMessage(
2075 "The object being loaded into does not match the triangulation "
2076 "that has been stored previously."));
2078 policy_name == ::internal::policy_to_string(*this->policy),
2079 ExcMessage("The policy currently associated with this DoFHandler (" +
2080 ::internal::policy_to_string(*this->policy) +
2081 ") does not match the one that was associated with the "
2082 "DoFHandler previously stored (" +
2083 policy_name + ")."));
2084 }
2085 else
2086 {
2087 ar &this->block_info_object;
2088 ar &number_cache;
2089
2090 object_dof_indices.clear();
2091
2092 object_dof_ptr.clear();
2093
2094 ar &this->object_dof_indices;
2095 ar &this->object_dof_ptr;
2096
2097 // these are the checks that correspond to the last block in the save()
2098 // function
2099 unsigned int n_cells;
2100 std::string fe_name;
2101 std::string policy_name;
2102
2103 ar &n_cells &fe_name &policy_name;
2104
2106 n_cells == this->tria->n_cells(),
2107 ExcMessage(
2108 "The object being loaded into does not match the triangulation "
2109 "that has been stored previously."));
2111 fe_name == this->get_fe(0).get_name(),
2112 ExcMessage(
2113 "The finite element associated with this DoFHandler does not match "
2114 "the one that was associated with the DoFHandler previously stored."));
2115 AssertThrow(policy_name == internal::policy_to_string(*this->policy),
2116 ExcMessage(
2117 "The policy currently associated with this DoFHandler (" +
2118 internal::policy_to_string(*this->policy) +
2119 ") does not match the one that was associated with the "
2120 "DoFHandler previously stored (" +
2121 policy_name + ")."));
2122 }
2123}
2124
2125
2126
2127template <int dim, int spacedim>
2131 const unsigned int level,
2132 const unsigned int dof_number,
2133 const unsigned int dofs_per_vertex)
2134{
2137 return indices[dofs_per_vertex * (level - coarsest_level) + dof_number];
2138}
2139
2140
2141
2142extern template class DoFHandler<1, 1>;
2143extern template class DoFHandler<1, 2>;
2144extern template class DoFHandler<1, 3>;
2145extern template class DoFHandler<2, 2>;
2146extern template class DoFHandler<2, 3>;
2147extern template class DoFHandler<3, 3>;
2148
2149
2150#endif // DOXYGEN
2151
2153
2154#endif
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition block_info.h:95
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:205
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:175
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
unsigned int level
Definition grid_out.cc:4626
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)
Definition exceptions.h:466
#define Assert(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:534
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:489
static ::ExceptionBase & ExcNoDominatedFiniteElementOnChildren()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:511
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:14871
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:59
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