Reference documentation for deal.II version 9.5.1
\(\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// Copyright (C) 1998 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_dof_handler_h
17#define dealii_dof_handler_h
18
19
20
21#include <deal.II/base/config.h>
22
28#include <deal.II/base/types.h>
29
36
38
39#include <boost/serialization/split_member.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#endif
98
315template <int dim, int spacedim = dim>
318{
323
324public:
337 using cell_accessor = typename ActiveSelector::CellAccessor;
338
351 using face_accessor = typename ActiveSelector::FaceAccessor;
352
360 using line_iterator = typename ActiveSelector::line_iterator;
361
375 using active_line_iterator = typename ActiveSelector::active_line_iterator;
376
384 using quad_iterator = typename ActiveSelector::quad_iterator;
385
399 using active_quad_iterator = typename ActiveSelector::active_quad_iterator;
400
408 using hex_iterator = typename ActiveSelector::hex_iterator;
409
419 using active_hex_iterator = typename ActiveSelector::active_hex_iterator;
420
441 using active_cell_iterator = typename ActiveSelector::active_cell_iterator;
442
469 using cell_iterator = typename ActiveSelector::cell_iterator;
470
487 using face_iterator = typename ActiveSelector::face_iterator;
488
500 using active_face_iterator = typename ActiveSelector::active_face_iterator;
501
502 using level_cell_accessor = typename LevelSelector::CellAccessor;
503 using level_face_accessor = typename LevelSelector::FaceAccessor;
504
505 using level_cell_iterator = typename LevelSelector::cell_iterator;
506 using level_face_iterator = typename LevelSelector::face_iterator;
507
508
512 static constexpr unsigned int dimension = dim;
513
517 static constexpr unsigned int space_dimension = spacedim;
518
522 static const types::fe_index default_fe_index = 0;
523
529 static const unsigned int invalid_fe_index DEAL_II_DEPRECATED =
531
538
542 using offset_type = unsigned int;
543
550 static const types::fe_index invalid_active_fe_index DEAL_II_DEPRECATED =
552
557 DoFHandler();
558
563
570 DoFHandler(const DoFHandler &) = delete;
571
575 virtual ~DoFHandler() override;
576
583 DoFHandler &
584 operator=(const DoFHandler &) = delete;
585
600 void
601 set_active_fe_indices(const std::vector<types::fe_index> &active_fe_indices);
602
608 DEAL_II_DEPRECATED_EARLY
609 void
610 set_active_fe_indices(const std::vector<unsigned int> &active_fe_indices);
611
625 std::vector<types::fe_index>
626 get_active_fe_indices() const;
627
644 DEAL_II_DEPRECATED_EARLY
645 void
646 get_active_fe_indices(std::vector<unsigned int> &active_fe_indices) const;
647
661 void
662 set_future_fe_indices(const std::vector<types::fe_index> &future_fe_indices);
663
674 std::vector<types::fe_index>
675 get_future_fe_indices() const;
676
684 void
685 reinit(const Triangulation<dim, spacedim> &tria);
686
716 void
717 distribute_dofs(const FiniteElement<dim, spacedim> &fe);
718
722 void
723 distribute_dofs(const hp::FECollection<dim, spacedim> &fe);
724
730 void
731 distribute_mg_dofs();
732
736 bool
738
744 bool
746
761 bool
763
770 void
771 initialize_local_block_info();
772
776 void
777 clear();
778
831 void
832 renumber_dofs(const std::vector<types::global_dof_index> &new_numbers);
833
838 void
839 renumber_dofs(const unsigned int level,
840 const std::vector<types::global_dof_index> &new_numbers);
841
870 unsigned int
871 max_couplings_between_dofs() const;
872
885 unsigned int
886 max_couplings_between_boundary_dofs() const;
887
888 /*--------------------------------------*/
889
894 /*
895 * @{
896 */
897
902 begin(const unsigned int level = 0) const;
903
921 begin_active(const unsigned int level = 0) const;
922
928 end() const;
929
935 end(const unsigned int level) const;
936
943 end_active(const unsigned int level) const;
944
950 begin_mg(const unsigned int level = 0) const;
951
957 end_mg(const unsigned int level) const;
958
964 end_mg() const;
965
981 cell_iterators() const;
982
1024 active_cell_iterators() const;
1025
1038 mg_cell_iterators() const;
1039
1056 cell_iterators_on_level(const unsigned int level) const;
1057
1074 active_cell_iterators_on_level(const unsigned int level) const;
1075
1092 mg_cell_iterators_on_level(const unsigned int level) const;
1093
1094 /*
1095 * @}
1096 */
1097
1098
1099 /*---------------------------------------*/
1100
1101
1126 n_dofs() const;
1127
1136 n_dofs(const unsigned int level) const;
1137
1143 n_boundary_dofs() const;
1144
1156 template <typename number>
1159 const std::map<types::boundary_id, const Function<spacedim, number> *>
1160 &boundary_ids) const;
1161
1167 n_boundary_dofs(const std::set<types::boundary_id> &boundary_ids) const;
1168
1184 const BlockInfo &
1185 block_info() const;
1186
1208
1214 const IndexSet &
1216
1221 const IndexSet &
1222 locally_owned_mg_dofs(const unsigned int level) const;
1223
1229 get_fe(const types::fe_index index = 0) const;
1230
1237
1243
1247 MPI_Comm
1249
1266 void
1267 prepare_for_serialization_of_active_fe_indices();
1268
1284 void
1285 deserialize_active_fe_indices();
1286
1295 virtual std::size_t
1296 memory_consumption() const;
1297
1303 template <class Archive>
1304 void
1305 save(Archive &ar, const unsigned int version) const;
1306
1312 template <class Archive>
1313 void
1314 load(Archive &ar, const unsigned int version);
1315
1316#ifdef DOXYGEN
1322 template <class Archive>
1323 void
1324 serialize(Archive &archive, const unsigned int version);
1325#else
1326 // This macro defines the serialize() method that is compatible with
1327 // the templated save() and load() method that have been implemented.
1328 BOOST_SERIALIZATION_SPLIT_MEMBER()
1329#endif
1330
1334 DeclException0(ExcNoFESelected);
1339 DeclException0(ExcInvalidBoundaryIndicator);
1344 DeclException1(ExcInvalidLevel,
1345 int,
1346 << "The given level " << arg1
1347 << " is not in the valid range!");
1352 DeclException1(ExcNewNumbersNotConsecutive,
1354 << "The given list of new dof indices is not consecutive: "
1355 << "the index " << arg1 << " does not exist.");
1359 DeclException2(ExcInvalidFEIndex,
1360 int,
1361 int,
1362 << "The mesh contains a cell with an active FE index of "
1363 << arg1 << ", but the finite element collection only has "
1364 << arg2 << " elements");
1365
1370 DeclExceptionMsg(ExcOnlyAvailableWithHP,
1371 "The current function doesn't make sense when used with a "
1372 "DoFHandler without hp-capabilities.");
1373
1378 DeclExceptionMsg(ExcNotImplementedWithHP,
1379 "The current function has not yet been implemented for a "
1380 "DoFHandler with hp-capabilities.");
1381
1382private:
1390 {
1391 public:
1395 MGVertexDoFs();
1396
1402 void
1403 init(const unsigned int coarsest_level,
1404 const unsigned int finest_level,
1405 const unsigned int dofs_per_vertex);
1406
1410 unsigned int
1411 get_coarsest_level() const;
1412
1416 unsigned int
1417 get_finest_level() const;
1418
1424 access_index(const unsigned int level,
1425 const unsigned int dof_number,
1426 const unsigned int dofs_per_vertex);
1427
1428 private:
1432 unsigned int coarsest_level;
1433
1437 unsigned int finest_level;
1438
1448 std::unique_ptr<types::global_dof_index[]> indices;
1449 };
1450
1458 {
1463 std::map<const cell_iterator, const types::fe_index>
1465
1470 std::map<const cell_iterator, const types::fe_index> refined_cells_fe_index;
1471
1476 std::map<const cell_iterator, const types::fe_index>
1478
1484 std::vector<types::fe_index> active_fe_indices;
1485
1491 std::unique_ptr<
1492 parallel::distributed::
1493 CellDataTransfer<dim, spacedim, std::vector<types::fe_index>>>
1495 };
1496
1501
1507
1513
1520
1525 std::unique_ptr<::internal::DoFHandlerImplementation::Policy::
1526 PolicyBase<dim, spacedim>>
1528
1537
1541 std::vector<::internal::DoFHandlerImplementation::NumberCache>
1543
1549 mutable std::vector<std::array<std::vector<types::global_dof_index>, dim + 1>>
1551
1560 mutable std::vector<std::array<std::vector<offset_type>, dim + 1>>
1562
1568 mutable std::array<std::vector<types::fe_index>, dim + 1>
1570
1574 mutable std::array<std::vector<offset_type>, dim + 1> hp_object_fe_ptr;
1575
1580 mutable std::vector<std::vector<types::fe_index>> hp_cell_active_fe_indices;
1581
1586 mutable std::vector<std::vector<types::fe_index>> hp_cell_future_fe_indices;
1587
1592 std::vector<MGVertexDoFs> mg_vertex_dofs;
1593
1597 std::vector<
1598 std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel<dim>>>
1600
1604 std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces<dim>>
1606
1611 std::unique_ptr<ActiveFEIndexTransfer> active_fe_index_transfer;
1612
1617 std::vector<boost::signals2::connection> tria_listeners;
1618
1624 std::vector<boost::signals2::connection> tria_listeners_for_transfer;
1625
1629 void
1630 clear_space();
1631
1635 void
1636 clear_mg_space();
1637
1641 void
1642 setup_policy();
1643
1647 void
1648 connect_to_triangulation_signals();
1649
1661 void
1662 create_active_fe_table();
1663
1676 void
1677 update_active_fe_table();
1678
1688 void
1689 pre_transfer_action();
1690
1699 void
1700 post_transfer_action();
1701
1710 void
1711 pre_distributed_transfer_action();
1712
1721 void
1722 post_distributed_transfer_action();
1723
1724
1725 // Make accessor objects friends.
1726 template <int, int, int, bool>
1727 friend class ::DoFAccessor;
1728 template <int, int, bool>
1729 friend class ::DoFCellAccessor;
1730 friend struct ::internal::DoFAccessorImplementation::Implementation;
1731 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1732
1733 // Likewise for DoFLevel objects since they need to access the vertex dofs
1734 // in the functions that set and retrieve vertex dof indices.
1735 friend struct ::internal::DoFHandlerImplementation::Implementation;
1736 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1737 friend struct ::internal::DoFHandlerImplementation::Policy::
1738 Implementation;
1739
1740 // explicitly check for sensible template arguments, but not on windows
1741 // because MSVC creates bogus warnings during normal compilation
1742#ifndef DEAL_II_MSVC
1743 static_assert(dim <= spacedim,
1744 "The dimension <dim> of a DoFHandler must be less than or "
1745 "equal to the space dimension <spacedim> in which it lives.");
1746#endif
1747};
1748
1749namespace internal
1750{
1751 namespace hp
1752 {
1753 namespace DoFHandlerImplementation
1754 {
1766 template <int dim, int spacedim>
1767 void
1769
1794 template <int dim, int spacedim = dim>
1795 unsigned int
1797 const typename DoFHandler<dim, spacedim>::cell_iterator &parent);
1798
1804 "No FiniteElement has been found in your FECollection that is "
1805 "dominated by all children of a cell you are trying to coarsen!");
1806 } // namespace DoFHandlerImplementation
1807 } // namespace hp
1808} // namespace internal
1809
1810#ifndef DOXYGEN
1811
1812/* ----------------------- Inline functions ----------------------------------
1813 */
1814
1815
1816template <int dim, int spacedim>
1819{
1820 return hp_capability_enabled;
1821}
1822
1823
1824
1825template <int dim, int spacedim>
1828{
1829 return mg_number_cache.size() > 0;
1830}
1831
1832
1833
1834template <int dim, int spacedim>
1837{
1838 return number_cache.n_global_dofs > 0;
1839}
1840
1841
1842
1843template <int dim, int spacedim>
1846{
1847 return number_cache.n_global_dofs;
1848}
1849
1850
1851
1852template <int dim, int spacedim>
1855 DoFHandler<dim, spacedim>::n_dofs(const unsigned int level) const
1856{
1857 Assert(has_level_dofs(),
1858 ExcMessage(
1859 "n_dofs(level) can only be called after distribute_mg_dofs()"));
1860 Assert(level < mg_number_cache.size(), ExcInvalidLevel(level));
1861 return mg_number_cache[level].n_global_dofs;
1862}
1863
1864
1865
1866template <int dim, int spacedim>
1869{
1870 return number_cache.n_locally_owned_dofs;
1871}
1872
1873
1874
1875template <int dim, int spacedim>
1878{
1879 return number_cache.locally_owned_dofs;
1880}
1881
1882
1883
1884template <int dim, int spacedim>
1887 const unsigned int level) const
1888{
1889 Assert(level < this->get_triangulation().n_global_levels(),
1890 ExcMessage("The given level index exceeds the number of levels "
1891 "present in the triangulation"));
1892 Assert(
1893 mg_number_cache.size() == this->get_triangulation().n_global_levels(),
1894 ExcMessage(
1895 "The level dofs are not set up properly! Did you call distribute_mg_dofs()?"));
1896 return mg_number_cache[level].locally_owned_dofs;
1897}
1898
1899
1900
1901template <int dim, int spacedim>
1904 const types::fe_index number) const
1905{
1906 Assert(fe_collection.size() > 0,
1907 ExcMessage("No finite element collection is associated with "
1908 "this DoFHandler"));
1909 return fe_collection[number];
1910}
1911
1912
1913
1914template <int dim, int spacedim>
1918{
1919 return fe_collection;
1920}
1921
1922
1923
1924template <int dim, int spacedim>
1928{
1929 Assert(tria != nullptr,
1930 ExcMessage("This DoFHandler object has not been associated "
1931 "with a triangulation."));
1932 return *tria;
1933}
1934
1935
1936
1937template <int dim, int spacedim>
1940{
1941 Assert(tria != nullptr,
1942 ExcMessage("This DoFHandler object has not been associated "
1943 "with a triangulation."));
1944 return tria->get_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:96
types::global_dof_index & access_index(const unsigned int level, const unsigned int dof_number, const unsigned int dofs_per_vertex)
std::unique_ptr< types::global_dof_index[]> indices
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
std::vector< std::vector< types::fe_index > > hp_cell_future_fe_indices
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
const hp::FECollection< dim, spacedim > & get_fe_collection() const
void save(Archive &ar, const unsigned int version) const
BlockInfo block_info_object
bool has_active_dofs() const
std::vector< boost::signals2::connection > tria_listeners_for_transfer
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
types::global_dof_index n_boundary_dofs() const
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
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
const IndexSet & locally_owned_dofs() const
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)
types::fe_index active_fe_index_type
DoFHandler & operator=(const DoFHandler &)=delete
DoFHandler(const DoFHandler &)=delete
bool hp_capability_enabled
bool has_hp_capabilities() const
types::global_dof_index n_dofs() 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
std::array< std::vector< offset_type >, dim+1 > hp_object_fe_ptr
MPI_Comm get_communicator() const
types::global_dof_index n_boundary_dofs(const std::map< types::boundary_id, const Function< spacedim, number > * > &boundary_ids) const
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
types::global_dof_index n_locally_owned_dofs() const
const BlockInfo & block_info() const
::internal::DoFHandlerImplementation::NumberCache number_cache
virtual MPI_Comm get_communicator() const
Definition tria.cc:11241
unsigned int n_cells() const
Definition tria.cc:13875
#define DEAL_II_DEPRECATED
Definition config.h:172
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:160
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
unsigned int level
Definition grid_out.cc:4618
#define DeclException0(Exception0)
Definition exceptions.h:465
#define Assert(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:533
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
static ::ExceptionBase & ExcNoDominatedFiniteElementOnChildren()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
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:118
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition tria.cc:13826
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)
const types::fe_index invalid_fe_index
Definition types.h:228
unsigned short int fe_index
Definition types.h:60
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
const ::Triangulation< dim, spacedim > & tria