Reference documentation for deal.II version GIT c9976103bc 2022-12-09 17:30:02+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\}}\)
dof_handler.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2022 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 
24 #include <deal.II/base/function.h>
25 #include <deal.II/base/index_set.h>
28 #include <deal.II/base/types.h>
29 
32 #include <deal.II/dofs/dof_faces.h>
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
50 template <int dim, int spacedim>
51 class FiniteElement;
52 template <int dim, int spacedim>
53 class Triangulation;
54 
55 namespace 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 
88 namespace parallel
89 {
90  namespace distributed
91  {
92  template <int dim, int spacedim, typename VectorType>
93  class CellDataTransfer;
94  }
95 } // namespace parallel
96 #endif
97 
312 template <int dim, int spacedim = dim>
313 class DoFHandler : public Subscriptor
314 {
319 
320 public:
333  using cell_accessor = typename ActiveSelector::CellAccessor;
334 
347  using face_accessor = typename ActiveSelector::FaceAccessor;
348 
356  using line_iterator = typename ActiveSelector::line_iterator;
357 
371  using active_line_iterator = typename ActiveSelector::active_line_iterator;
372 
380  using quad_iterator = typename ActiveSelector::quad_iterator;
381 
395  using active_quad_iterator = typename ActiveSelector::active_quad_iterator;
396 
404  using hex_iterator = typename ActiveSelector::hex_iterator;
405 
415  using active_hex_iterator = typename ActiveSelector::active_hex_iterator;
416 
437  using active_cell_iterator = typename ActiveSelector::active_cell_iterator;
438 
465  using cell_iterator = typename ActiveSelector::cell_iterator;
466 
483  using face_iterator = typename ActiveSelector::face_iterator;
484 
496  using active_face_iterator = typename ActiveSelector::active_face_iterator;
497 
498  using level_cell_accessor = typename LevelSelector::CellAccessor;
499  using level_face_accessor = typename LevelSelector::FaceAccessor;
500 
501  using level_cell_iterator = typename LevelSelector::cell_iterator;
502  using level_face_iterator = typename LevelSelector::face_iterator;
503 
504 
508  static constexpr unsigned int dimension = dim;
509 
513  static constexpr unsigned int space_dimension = spacedim;
514 
519 
525  static const unsigned int invalid_fe_index DEAL_II_DEPRECATED =
527 
534 
538  using offset_type = unsigned int;
539 
548 
554 
559 
566  DoFHandler(const DoFHandler &) = delete;
567 
571  virtual ~DoFHandler() override;
572 
579  DoFHandler &
580  operator=(const DoFHandler &) = delete;
581 
595  void
596  set_active_fe_indices(const std::vector<types::fe_index> &active_fe_indices);
597 
603  DEAL_II_DEPRECATED_EARLY
604  void
605  set_active_fe_indices(const std::vector<unsigned int> &active_fe_indices);
606 
620  std::vector<types::fe_index>
622 
639  DEAL_II_DEPRECATED_EARLY
640  void
641  get_active_fe_indices(std::vector<unsigned int> &active_fe_indices) const;
642 
656  void
657  set_future_fe_indices(const std::vector<types::fe_index> &future_fe_indices);
658 
669  std::vector<types::fe_index>
671 
679  void
681 
711  void
713 
717  void
719 
725  void
727 
731  bool
733 
739  bool
740  has_level_dofs() const;
741 
756  bool
758 
765  void
767 
771  void
772  clear();
773 
826  void
827  renumber_dofs(const std::vector<types::global_dof_index> &new_numbers);
828 
833  void
834  renumber_dofs(const unsigned int level,
835  const std::vector<types::global_dof_index> &new_numbers);
836 
865  unsigned int
867 
880  unsigned int
882 
883  /*--------------------------------------*/
884 
889  /*
890  * @{
891  */
892 
897  begin(const unsigned int level = 0) const;
898 
916  begin_active(const unsigned int level = 0) const;
917 
923  end() const;
924 
930  end(const unsigned int level) const;
931 
938  end_active(const unsigned int level) const;
939 
945  begin_mg(const unsigned int level = 0) const;
946 
952  end_mg(const unsigned int level) const;
953 
959  end_mg() const;
960 
976  cell_iterators() const;
977 
1020 
1034 
1051  cell_iterators_on_level(const unsigned int level) const;
1052 
1069  active_cell_iterators_on_level(const unsigned int level) const;
1070 
1087  mg_cell_iterators_on_level(const unsigned int level) const;
1088 
1089  /*
1090  * @}
1091  */
1092 
1093 
1094  /*---------------------------------------*/
1095 
1096 
1121  n_dofs() const;
1122 
1131  n_dofs(const unsigned int level) const;
1132 
1139 
1151  template <typename number>
1154  const std::map<types::boundary_id, const Function<spacedim, number> *>
1155  &boundary_ids) const;
1156 
1162  n_boundary_dofs(const std::set<types::boundary_id> &boundary_ids) const;
1163 
1179  const BlockInfo &
1180  block_info() const;
1181 
1203 
1209  const IndexSet &
1211 
1216  const IndexSet &
1217  locally_owned_mg_dofs(const unsigned int level) const;
1218 
1224  get_fe(const types::fe_index index = 0) const;
1225 
1232 
1238 
1242  MPI_Comm
1244 
1261  void
1263 
1279  void
1281 
1290  virtual std::size_t
1292 
1298  template <class Archive>
1299  void
1300  save(Archive &ar, const unsigned int version) const;
1301 
1307  template <class Archive>
1308  void
1309  load(Archive &ar, const unsigned int version);
1310 
1311 #ifdef DOXYGEN
1317  template <class Archive>
1318  void
1319  serialize(Archive &archive, const unsigned int version);
1320 #else
1321  // This macro defines the serialize() method that is compatible with
1322  // the templated save() and load() method that have been implemented.
1323  BOOST_SERIALIZATION_SPLIT_MEMBER()
1324 #endif
1325 
1340  int,
1341  << "The given level " << arg1
1342  << " is not in the valid range!");
1349  << "The given list of new dof indices is not consecutive: "
1350  << "the index " << arg1 << " does not exist.");
1355  int,
1356  int,
1357  << "The mesh contains a cell with an active FE index of "
1358  << arg1 << ", but the finite element collection only has "
1359  << arg2 << " elements");
1360 
1366  "The current function doesn't make sense when used with a "
1367  "DoFHandler without hp-capabilities.");
1368 
1374  "The current function has not yet been implemented for a "
1375  "DoFHandler with hp-capabilities.");
1376 
1377 private:
1385  {
1386  public:
1391 
1397  void
1398  init(const unsigned int coarsest_level,
1399  const unsigned int finest_level,
1400  const unsigned int dofs_per_vertex);
1401 
1405  unsigned int
1407 
1411  unsigned int
1413 
1419  access_index(const unsigned int level,
1420  const unsigned int dof_number,
1421  const unsigned int dofs_per_vertex);
1422 
1423  private:
1427  unsigned int coarsest_level;
1428 
1432  unsigned int finest_level;
1433 
1443  std::unique_ptr<types::global_dof_index[]> indices;
1444  };
1445 
1453  {
1458  std::map<const cell_iterator, const types::fe_index>
1460 
1465  std::map<const cell_iterator, const types::fe_index> refined_cells_fe_index;
1466 
1471  std::map<const cell_iterator, const types::fe_index>
1473 
1479  std::vector<types::fe_index> active_fe_indices;
1480 
1486  std::unique_ptr<
1487  parallel::distributed::
1488  CellDataTransfer<dim, spacedim, std::vector<types::fe_index>>>
1490  };
1491 
1496 
1502 
1508 
1515 
1520  std::unique_ptr<::internal::DoFHandlerImplementation::Policy::
1521  PolicyBase<dim, spacedim>>
1523 
1532 
1536  std::vector<::internal::DoFHandlerImplementation::NumberCache>
1538 
1544  mutable std::vector<std::array<std::vector<types::global_dof_index>, dim + 1>>
1546 
1555  mutable std::vector<std::array<std::vector<offset_type>, dim + 1>>
1557 
1563  mutable std::array<std::vector<types::fe_index>, dim + 1>
1565 
1569  mutable std::array<std::vector<offset_type>, dim + 1> hp_object_fe_ptr;
1570 
1575  mutable std::vector<std::vector<types::fe_index>> hp_cell_active_fe_indices;
1576 
1581  mutable std::vector<std::vector<types::fe_index>> hp_cell_future_fe_indices;
1582 
1587  std::vector<MGVertexDoFs> mg_vertex_dofs;
1588 
1592  std::vector<
1593  std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel<dim>>>
1595 
1599  std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces<dim>>
1601 
1606  std::unique_ptr<ActiveFEIndexTransfer> active_fe_index_transfer;
1607 
1612  std::vector<boost::signals2::connection> tria_listeners;
1613 
1619  std::vector<boost::signals2::connection> tria_listeners_for_transfer;
1620 
1624  void
1626 
1630  void
1632 
1636  void
1638 
1642  void
1644 
1656  void
1658 
1671  void
1673 
1683  void
1685 
1694  void
1696 
1705  void
1707 
1716  void
1718 
1719 
1720  // Make accessor objects friends.
1721  template <int, int, int, bool>
1722  friend class ::DoFAccessor;
1723  template <int, int, bool>
1724  friend class ::DoFCellAccessor;
1725  friend struct ::internal::DoFAccessorImplementation::Implementation;
1726  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1727 
1728  // Likewise for DoFLevel objects since they need to access the vertex dofs
1729  // in the functions that set and retrieve vertex dof indices.
1730  friend struct ::internal::DoFHandlerImplementation::Implementation;
1731  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1732  friend struct ::internal::DoFHandlerImplementation::Policy::
1733  Implementation;
1734 
1735  // explicitly check for sensible template arguments, but not on windows
1736  // because MSVC creates bogus warnings during normal compilation
1737 #ifndef DEAL_II_MSVC
1738  static_assert(dim <= spacedim,
1739  "The dimension <dim> of a DoFHandler must be less than or "
1740  "equal to the space dimension <spacedim> in which it lives.");
1741 #endif
1742 };
1743 
1744 namespace internal
1745 {
1746  namespace hp
1747  {
1748  namespace DoFHandlerImplementation
1749  {
1761  template <int dim, int spacedim>
1762  void
1764 
1789  template <int dim, int spacedim = dim>
1790  unsigned int
1792  const typename DoFHandler<dim, spacedim>::cell_iterator &parent);
1793 
1799  "No FiniteElement has been found in your FECollection that is "
1800  "dominated by all children of a cell you are trying to coarsen!");
1801  } // namespace DoFHandlerImplementation
1802  } // namespace hp
1803 } // namespace internal
1804 
1805 #ifndef DOXYGEN
1806 
1807 /* ----------------------- Inline functions ----------------------------------
1808  */
1809 
1810 
1811 template <int dim, int spacedim>
1812 inline bool
1814 {
1815  return hp_capability_enabled;
1816 }
1817 
1818 
1819 
1820 template <int dim, int spacedim>
1821 inline bool
1823 {
1824  return mg_number_cache.size() > 0;
1825 }
1826 
1827 
1828 
1829 template <int dim, int spacedim>
1830 inline bool
1832 {
1833  return number_cache.n_global_dofs > 0;
1834 }
1835 
1836 
1837 
1838 template <int dim, int spacedim>
1841 {
1842  return number_cache.n_global_dofs;
1843 }
1844 
1845 
1846 
1847 template <int dim, int spacedim>
1849 DoFHandler<dim, spacedim>::n_dofs(const unsigned int level) const
1850 {
1852  ExcMessage(
1853  "n_dofs(level) can only be called after distribute_mg_dofs()"));
1855  return mg_number_cache[level].n_global_dofs;
1856 }
1857 
1858 
1859 
1860 template <int dim, int spacedim>
1863 {
1865 }
1866 
1867 
1868 
1869 template <int dim, int spacedim>
1870 const IndexSet &
1872 {
1874 }
1875 
1876 
1877 
1878 template <int dim, int spacedim>
1879 const IndexSet &
1881 {
1882  Assert(level < this->get_triangulation().n_global_levels(),
1883  ExcMessage("The given level index exceeds the number of levels "
1884  "present in the triangulation"));
1885  Assert(
1886  mg_number_cache.size() == this->get_triangulation().n_global_levels(),
1887  ExcMessage(
1888  "The level dofs are not set up properly! Did you call distribute_mg_dofs()?"));
1889  return mg_number_cache[level].locally_owned_dofs;
1890 }
1891 
1892 
1893 
1894 template <int dim, int spacedim>
1895 inline const FiniteElement<dim, spacedim> &
1897 {
1898  Assert(fe_collection.size() > 0,
1899  ExcMessage("No finite element collection is associated with "
1900  "this DoFHandler"));
1901  return fe_collection[number];
1902 }
1903 
1904 
1905 
1906 template <int dim, int spacedim>
1907 inline const hp::FECollection<dim, spacedim> &
1909 {
1910  return fe_collection;
1911 }
1912 
1913 
1914 
1915 template <int dim, int spacedim>
1916 inline const Triangulation<dim, spacedim> &
1918 {
1919  Assert(tria != nullptr,
1920  ExcMessage("This DoFHandler object has not been associated "
1921  "with a triangulation."));
1922  return *tria;
1923 }
1924 
1925 
1926 
1927 template <int dim, int spacedim>
1928 inline MPI_Comm
1930 {
1931  Assert(tria != nullptr,
1932  ExcMessage("This DoFHandler object has not been associated "
1933  "with a triangulation."));
1934  return tria->get_communicator();
1935 }
1936 
1937 
1938 
1939 template <int dim, int spacedim>
1940 inline const BlockInfo &
1942 {
1944 
1945  return block_info_object;
1946 }
1947 
1948 
1949 
1950 template <int dim, int spacedim>
1951 template <typename number>
1954  const std::map<types::boundary_id, const Function<spacedim, number> *>
1955  &boundary_ids) const
1956 {
1957  Assert(!(dim == 2 && spacedim == 3) || this->hp_capability_enabled == false,
1959 
1960  // extract the set of boundary ids and forget about the function object
1961  // pointers
1962  std::set<types::boundary_id> boundary_ids_only;
1963  for (typename std::map<types::boundary_id,
1964  const Function<spacedim, number> *>::const_iterator p =
1965  boundary_ids.begin();
1966  p != boundary_ids.end();
1967  ++p)
1968  boundary_ids_only.insert(p->first);
1969 
1970  // then just hand everything over to the other function that does the work
1971  return n_boundary_dofs(boundary_ids_only);
1972 }
1973 
1974 
1975 
1976 namespace internal
1977 {
1985  template <int dim, int spacedim>
1986  std::string
1987  policy_to_string(const ::internal::DoFHandlerImplementation::Policy::
1988  PolicyBase<dim, spacedim> &policy);
1989 } // namespace internal
1990 
1991 
1992 
1993 template <int dim, int spacedim>
1994 template <class Archive>
1995 void
1996 DoFHandler<dim, spacedim>::save(Archive &ar, const unsigned int) const
1997 {
1998  if (this->hp_capability_enabled)
1999  {
2000  ar &this->object_dof_indices;
2001  ar &this->object_dof_ptr;
2002 
2003  ar &this->hp_cell_active_fe_indices;
2004  ar &this->hp_cell_future_fe_indices;
2005 
2006  ar &hp_object_fe_ptr;
2008 
2009  ar &number_cache;
2010 
2011  ar &mg_number_cache;
2012 
2013  // write out the number of triangulation cells and later check during
2014  // loading that this number is indeed correct; same with something that
2015  // identifies the policy
2016  const unsigned int n_cells = this->tria->n_cells();
2017  std::string policy_name =
2019 
2020  ar &n_cells &policy_name;
2021  }
2022  else
2023  {
2024  ar &this->block_info_object;
2025  ar &number_cache;
2026 
2027  ar &this->object_dof_indices;
2028  ar &this->object_dof_ptr;
2029 
2030  // write out the number of triangulation cells and later check during
2031  // loading that this number is indeed correct; same with something that
2032  // identifies the FE and the policy
2033  unsigned int n_cells = this->tria->n_cells();
2034  std::string fe_name = this->get_fe(0).get_name();
2035  std::string policy_name = internal::policy_to_string(*this->policy);
2036 
2037  ar &n_cells &fe_name &policy_name;
2038  }
2039 }
2040 
2041 
2042 
2043 template <int dim, int spacedim>
2044 template <class Archive>
2045 void
2046 DoFHandler<dim, spacedim>::load(Archive &ar, const unsigned int)
2047 {
2048  if (this->hp_capability_enabled)
2049  {
2050  ar &this->object_dof_indices;
2051  ar &this->object_dof_ptr;
2052 
2053  ar &this->hp_cell_active_fe_indices;
2054  ar &this->hp_cell_future_fe_indices;
2055 
2056  ar &hp_object_fe_ptr;
2058 
2059  ar &number_cache;
2060 
2061  ar &mg_number_cache;
2062 
2063  // these are the checks that correspond to the last block in the save()
2064  // function
2065  unsigned int n_cells;
2066  std::string policy_name;
2067 
2068  ar &n_cells &policy_name;
2069 
2070  AssertThrow(
2071  n_cells == this->tria->n_cells(),
2072  ExcMessage(
2073  "The object being loaded into does not match the triangulation "
2074  "that has been stored previously."));
2075  AssertThrow(
2076  policy_name == ::internal::policy_to_string(*this->policy),
2077  ExcMessage("The policy currently associated with this DoFHandler (" +
2079  ") does not match the one that was associated with the "
2080  "DoFHandler previously stored (" +
2081  policy_name + ")."));
2082  }
2083  else
2084  {
2085  ar &this->block_info_object;
2086  ar &number_cache;
2087 
2088  object_dof_indices.clear();
2089 
2090  object_dof_ptr.clear();
2091 
2092  ar &this->object_dof_indices;
2093  ar &this->object_dof_ptr;
2094 
2095  // these are the checks that correspond to the last block in the save()
2096  // function
2097  unsigned int n_cells;
2098  std::string fe_name;
2099  std::string policy_name;
2100 
2101  ar &n_cells &fe_name &policy_name;
2102 
2103  AssertThrow(
2104  n_cells == this->tria->n_cells(),
2105  ExcMessage(
2106  "The object being loaded into does not match the triangulation "
2107  "that has been stored previously."));
2108  AssertThrow(
2109  fe_name == this->get_fe(0).get_name(),
2110  ExcMessage(
2111  "The finite element associated with this DoFHandler does not match "
2112  "the one that was associated with the DoFHandler previously stored."));
2113  AssertThrow(policy_name == internal::policy_to_string(*this->policy),
2114  ExcMessage(
2115  "The policy currently associated with this DoFHandler (" +
2117  ") does not match the one that was associated with the "
2118  "DoFHandler previously stored (" +
2119  policy_name + ")."));
2120  }
2121 }
2122 
2123 
2124 
2125 template <int dim, int spacedim>
2126 inline types::global_dof_index &
2128  const unsigned int level,
2129  const unsigned int dof_number,
2130  const unsigned int dofs_per_vertex)
2131 {
2132  Assert((level >= coarsest_level) && (level <= finest_level),
2133  ExcInvalidLevel(level));
2134  return indices[dofs_per_vertex * (level - coarsest_level) + dof_number];
2135 }
2136 
2137 
2138 
2139 extern template class DoFHandler<1, 1>;
2140 extern template class DoFHandler<1, 2>;
2141 extern template class DoFHandler<1, 3>;
2142 extern template class DoFHandler<2, 2>;
2143 extern template class DoFHandler<2, 3>;
2144 extern template class DoFHandler<3, 3>;
2145 
2146 
2147 #endif // DOXYGEN
2148 
2150 
2151 #endif
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:94
types::global_dof_index & access_index(const unsigned int level, const unsigned int dof_number, const unsigned int dofs_per_vertex)
unsigned int coarsest_level
Definition: dof_handler.h:1427
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
Definition: dof_handler.h:1443
unsigned int get_coarsest_level() const
cell_iterator end() const
void pre_distributed_transfer_action()
void post_transfer_action()
virtual std::size_t memory_consumption() const
unsigned int max_couplings_between_dofs() const
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > mg_levels
Definition: dof_handler.h:1594
std::vector< types::fe_index > get_active_fe_indices() const
std::unique_ptr< ActiveFEIndexTransfer > active_fe_index_transfer
Definition: dof_handler.h:1606
hp::FECollection< dim, spacedim > fe_collection
Definition: dof_handler.h:1514
void create_active_fe_table()
static constexpr unsigned int dimension
Definition: dof_handler.h:508
static const types::fe_index default_fe_index
Definition: dof_handler.h:518
std::vector< std::vector< types::fe_index > > hp_cell_future_fe_indices
Definition: dof_handler.h:1581
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
Definition: dof_handler.h:1507
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
Definition: dof_handler.h:1495
bool has_active_dofs() const
std::vector< boost::signals2::connection > tria_listeners_for_transfer
Definition: dof_handler.h:1619
level_cell_iterator end_mg() const
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
static constexpr unsigned int space_dimension
Definition: dof_handler.h:513
std::vector< types::fe_index > get_future_fe_indices() 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
Definition: dof_handler.h:1575
unsigned int offset_type
Definition: dof_handler.h:538
std::vector< MGVertexDoFs > mg_vertex_dofs
Definition: dof_handler.h:1587
std::vector< boost::signals2::connection > tria_listeners
Definition: dof_handler.h:1612
void post_distributed_transfer_action()
std::vector< std::array< std::vector< types::global_dof_index >, dim+1 > > object_dof_indices
Definition: dof_handler.h:1545
typename LevelSelector::CellAccessor level_cell_accessor
Definition: dof_handler.h:498
const BlockInfo & block_info() const
const Triangulation< dim, spacedim > & get_triangulation() const
void pre_transfer_action()
void clear_space()
void reinit(const Triangulation< dim, spacedim > &tria)
std::vector< std::array< std::vector< offset_type >, dim+1 > > object_dof_ptr
Definition: dof_handler.h:1556
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > mg_faces
Definition: dof_handler.h:1600
void serialize(Archive &archive, const unsigned int version)
types::fe_index active_fe_index_type
Definition: dof_handler.h:533
DoFHandler(const DoFHandler &)=delete
active_cell_iterator begin_active(const unsigned int level=0) const
void distribute_mg_dofs()
active_cell_iterator end_active(const unsigned int level) const
bool hp_capability_enabled
Definition: dof_handler.h:1501
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
Definition: dof_handler.h:501
void prepare_for_serialization_of_active_fe_indices()
std::vector<::internal::DoFHandlerImplementation::NumberCache > mg_number_cache
Definition: dof_handler.h:1537
void load(Archive &ar, const unsigned int version)
DoFHandler & operator=(const DoFHandler &)=delete
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
typename LevelSelector::face_iterator level_face_iterator
Definition: dof_handler.h:502
cell_iterator begin(const unsigned int level=0) const
static const types::fe_index invalid_active_fe_index
Definition: dof_handler.h:546
std::array< std::vector< offset_type >, dim+1 > hp_object_fe_ptr
Definition: dof_handler.h:1569
MPI_Comm get_communicator() const
static const unsigned int invalid_fe_index
Definition: dof_handler.h:525
void clear()
const IndexSet & locally_owned_dofs() const
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
bool has_level_dofs() const
std::unique_ptr<::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > > policy
Definition: dof_handler.h:1522
typename LevelSelector::FaceAccessor level_face_accessor
Definition: dof_handler.h:499
types::global_dof_index n_dofs(const unsigned int level) const
std::array< std::vector< types::fe_index >, dim+1 > hp_object_fe_indices
Definition: dof_handler.h:1564
const hp::FECollection< dim, spacedim > & get_fe_collection() const
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
void deserialize_active_fe_indices()
types::global_dof_index n_boundary_dofs(const std::set< types::boundary_id > &boundary_ids) const
::internal::DoFHandlerImplementation::NumberCache number_cache
Definition: dof_handler.h:1531
void set_future_fe_indices(const std::vector< types::fe_index > &future_fe_indices)
virtual ~DoFHandler() override
virtual std::string get_name() const =0
virtual MPI_Comm get_communicator() const
unsigned int n_cells() const
unsigned int size() const
Definition: collection.h:264
#define DEAL_II_DEPRECATED
Definition: config.h:164
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
unsigned int level
Definition: grid_out.cc:4608
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
IteratorRange< active_cell_iterator > active_cell_iterators() const
IteratorRange< level_cell_iterator > mg_cell_iterators() const
IteratorRange< cell_iterator > cell_iterators() const
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcInvalidFEIndex(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidBoundaryIndicator()
static ::ExceptionBase & ExcInvalidLevel(int arg1)
static ::ExceptionBase & ExcNewNumbersNotConsecutive(types::global_dof_index arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1501
static ::ExceptionBase & ExcNoFESelected()
static ::ExceptionBase & ExcOnlyAvailableWithHP()
static ::ExceptionBase & ExcNotImplementedWithHP()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcNoDominatedFiniteElementOnChildren()
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1611
typename ActiveSelector::FaceAccessor face_accessor
Definition: dof_handler.h:347
typename ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:465
typename ActiveSelector::line_iterator line_iterator
Definition: dof_handler.h:356
typename ActiveSelector::face_iterator face_iterator
Definition: dof_handler.h:483
typename ActiveSelector::active_line_iterator active_line_iterator
Definition: dof_handler.h:371
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:437
typename ActiveSelector::CellAccessor cell_accessor
Definition: dof_handler.h:333
typename ActiveSelector::quad_iterator quad_iterator
Definition: dof_handler.h:380
typename ActiveSelector::active_hex_iterator active_hex_iterator
Definition: dof_handler.h:415
typename ActiveSelector::active_quad_iterator active_quad_iterator
Definition: dof_handler.h:395
typename ActiveSelector::hex_iterator hex_iterator
Definition: dof_handler.h:404
typename ActiveSelector::active_face_iterator active_face_iterator
Definition: dof_handler.h:496
Definition: hp.h:118
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:13750
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)
Definition: dof_handler.cc:53
const types::fe_index invalid_fe_index
Definition: types.h:221
unsigned int global_dof_index
Definition: types.h:81
unsigned short int fe_index
Definition: types.h:59
unsigned int boundary_id
Definition: types.h:134
std::map< const cell_iterator, const types::fe_index > refined_cells_fe_index
Definition: dof_handler.h:1465
std::map< const cell_iterator, const types::fe_index > persisting_cells_fe_index
Definition: dof_handler.h:1459
std::vector< types::fe_index > active_fe_indices
Definition: dof_handler.h:1479
std::map< const cell_iterator, const types::fe_index > coarsened_cells_fe_index
Definition: dof_handler.h:1472
std::unique_ptr< parallel::distributed::CellDataTransfer< dim, spacedim, std::vector< types::fe_index > > > cell_data_transfer
Definition: dof_handler.h:1489